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ABSTRACT 


We report the discovery of an accreting supermassive black hole at z = 8.679. This galaxy, denoted here 
as CEERS 1019, was previously discovered as a Lya-break galaxy by Hubble with a Lya redshift from Keck. 
As part of the Cosmic Evolution Early Release Science (CEERS) survey, we have observed this source with 
JWST/NIRSpec spectroscopy, MIRI imaging, NIRCam imaging, and NIRCam/WESS slitless spectroscopy. 
We use the NIRSpec R ~ 1000 medium-resolution grating spectra covering 1—5ym to uncover a plethora of 
significantly-detected emission lines. The extremely strong [O III] emission line, observed at 4.8473 jum, 
confirms the ground-based Lya redshift. We detect a significant broad (FWHM ~ 1200 km s^!) component 
in the H8 emission line, which we conclude originates in the broad-line region of an active galactic nucleus 
(AGN) as the lack of a similar broad component in other, stronger, forbidden lines rejects an outflow origin. This 
hypothesis is supported by the presence of weak high-ionization lines (N V, N IV], C IIIJ), some of which have 
broad components, as well as a spatial point-source component embedded within a smoother surface brightness 
profile. The implied mass of the black hole is log (Mga /Mo) = 6.95 + 0.37, and we estimate that it is accreting 
at 1.2 (+ 0.5) x the Eddington limit. The 1-8 jum photometric spectral energy distribution (SED) from NIRCam 
and MIRI shows a continuum dominated by starlight and constrains the host galaxy to be massive (log M/Ms ~ 
9.5) and highly star-forming (SFR ~ 30 Mo yr^!; log sSFR ~ —7.9 yr^!). The ratios of detected strong emission 
lines show that the gas in this galaxy is metal-poor (Z/Z 0.1), dense (n, ~ 10? cm^?), and highly ionized (log 
U ~ —2.1), consistent with the general galaxy population observed with JWST at such high redshifts that do not 
display significant AGN signatures. We use this presently highest-redshift AGN discovery to place constraints 
on black hole seeding models and find that a combination of either super-Eddington accretion from stellar seeds 
or Eddington accretion from very massive black hole seeds is required to form this object by the observed epoch. 


1. INTRODUCTION 


One of the most consequential periods in cosmic history is 
the Epoch of Reionization (EoR), where the material between 
galaxies underwent a significant transition from neutral to 
ionized hydrogen. Satisfactory explanations for the sources of 
radiation that contributed to this process, when it started, and 
how long it lasted are all presently wanting, primarily from 
the lack of necessary observatories and instruments tuned to 
the early Universe — until now. With the launch of JWST, 
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we are on the cusp of an exciting new era in astronomy as 
detailed studies of galaxies in the first billion years are finally 
possible. 

One of the key questions JWST was designed to answer 
was when and how the first black holes formed. Defining 
this epoch will help constrain the role these sources played in 
reionization alongside ionizing photons from massive stars. 
Supermassive black holes (SMBHs) exist at the centers of 
massive galaxies in the present-day Universe and exhibit a 
well-studied correlation with the velocity dispersions (and 
stellar masses) of galaxy bulge components (see review by 
Kormendy & Ho 2013). These massive objects have >13 
Gyr of cosmic time to grow to their present-day mass, which 
is possible via standard accretion scenarios with stellar mass 
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black hole seeds (~1-10 Mo; possibly up to 100 Mo if a 
Population III star remnant). However, the surprising dis- 
covery of z > 6 quasars with the Sloan Digital Sky Survey 
(SDSS) with black hole masses of log (MgxH/Mo) > 9 (e.g., 
Fan et al. 2006; Bafiados et al. 2018; Wang et al. 2021; Farina 
et al. 2022) challenges models of black hole growth. Such 
objects require very early (z ~ 25-30) stellar mass seeds with 
near-unity duty cycles of Eddington-limited accretion and/or 
super-Eddington accretion to grow to such a mass from a 
stellar mass seed in < 1 Gyr (e.g., Volonteri et al. 2021 and 
references therein). None of these scenarios are adequately 
predicted by simulations (e.g., Volonteri et al. 2021; Fontanot 
et al. 2023), and this tension is further exacerbated by sim- 
ilarly massive quasars now being discovered at z > 7 (e.g., 
Mortlock et al. 2011; Bafiados et al. 2018, see Fan et al. 2022 
for a recent review), with the current highest-redshift known 
AGN being a bright quasar at z = 7.64 (Wang et al. 2021). 

These surprisingly-massive black holes in the first Gyr of 
cosmic history have led to an alternative seeding theory: di- 
rect collapse black holes (DCBHs; Bromm & Loeb 2003). 
In this model, mini-halos irradiated by Lyman-Werner pho- 
tons (11.2-13.6 eV) cannot form molecular gas and thus do 
not form Population III stars. As these halos grow, they 
eventually cross the atomic cooling regime, at which point 
the gas will begin to cool rapidly via H I cooling (primarily 
Lya). Further fragmentation-inducing cooling is avoided due 
to the H2-suppressing UV background, leading to collapse 
into massive black holes in the range of 10^7 Mo, preceded 
by a brief phase as a supermassive star (see Smith & Bromm 
2019, Woods et al. 2019 and references therein for a detailed 
explanation of this process). A similar massive black hole 
seed could be formed as the remnant of a supermassive star 
powered by WIMP-like dark matter annihilation (so-called 
“dark stars"; Freese et al. 2016; Ilie et al. 2012). Such mas- 
sive black hole seeds, forming at z ~ 10-15 (by necessity 
after the first generation of stars), could more feasibly grow 
into the observed z ~ 6-7 quasar population (e.g., Madau 
et al. 2014; Natarajan et al. 2017; Regan et al. 2019; Latif 
et al. 2021; Pacucci & Loeb 2022; Trinca et al. 2023; Mas- 
sonneau et al. 2023). These massive quasars must represent 
only the extreme cases — there likely exists a much larger 
population of lower-mass black holes, and/or obscured black 
holes, waiting to be discovered. 

While DCBHs could alleviate tension between observed 
black hole masses and our theories of black hole growth, 
such objects have yet to be observed. One clear pathway to 
understanding black hole growth is to observe more SMBHs 
in the epoch of reionization. Identifying modest-sized black 
holes at earlier cosmic times could provide further evidence as 
to whether DCBHs are a necessary pathway. Additionally, the 
discovery of such a population would both better explain how 
the observed z ~ 6 quasar population originally formed and 


inform on the potential contribution of AGNS to reionization 
both through X-ray heating (e.g., Jeon et al. 2014) and through 
ionizing photon contributions (e.g., Finkelstein et al. 2019; 
Giallongo et al. 2019; Grazian et al. 2020; Yung et al. 2021; 
Grazian et al. 2022). 

Prior to JWST, only the most massive SMBHs at high red- 
shifts could be identified. However, the spectroscopic ca- 
pabilities of JWST now enable the search for signs of AGN 
activity from less luminous sources, particularly those em- 
bedded in galaxies whose stellar emission dominates the total 
galaxy luminosity (e.g., Endsley et al. 2022; Kocevski et al. 
2023), and/or where the bulk of the accretion emission is 
obscured (e.g., Fujimoto et al. 2022, see also Furtak et al. 
2022). 

Here we report the discovery of the first known AGN at 
z » 8. The galaxy harboring this AGN was first identi- 
fied as a candidate z ~ 8 galaxy by Roberts-Borsani et al. 
(2016a, named EGSY-2008532660). Its spectroscopic red- 
shift was measured via Lya emission via Keck/MOSFIRE to 
be zpyo = 8.683*000! from (Zitrin et al. 2015; as EGSY8p7); 
at the time, and for several years, was the farthest known 
Lya-emitter. It later gained a potential Nv detection at 
Zsys = 8.667 from Mainali et al. (2018), from additional 
Keck/MOSFIRE data, and was also the brightest z > 8.5 
galaxy found in the CANDELS survey (Finkelstein et al. 
2022a; EGS, 7910 6811). Here we report the results from 
JWST with the CEERS (Finkelstein et al. 2022b) dataset. The 
NIRSpec spectroscopy (Arrabal Haro et al., in prep) of this 
source has an ID from the micro-shutter array (MSA) of 1019; 
thus, and hereafter, we refer to it as CEERS. 1019. 

In 82, we present the data from four different JWST obser- 
vational modes from the CEERS Survey. In $3, we describe 
our method of emission line fitting while we explore the de- 
tected emission lines in detail in $4. In $5, we analyze the 
properties of this galaxy from the available imaging data and 
from our spectroscopic data in $6. We discuss the implica- 
tions of these results in $7 and present our conclusions in 
$8. Throughout this paper we assume a flat Planck cosmol- 
ogy with Ho = 67.36 km s7! Mpc^!, €, = 0.3153 and Qa = 
0.6847 (Planck Collaboration et al. 2020). All magnitudes are 
in the AB system, and all rest-frame wavelengths are vacuum. 


2. DATA 


Data presented in this work were taken as part of the Cosmic 
Evolution Early Release Science Survey (CEERS; ERS 1345, 
PI: S. Finkelstein, Finkelstein et al. 2022b; Bagley et al. 2022) 
inthe CANDELS (Grogin et al. 2011; Koekemoer et al. 2011) 
Extended Groth Strip (EGS) field. The complete details of 
the CEERS program will be presented in Finkelstein et al., 
(in prep.). This source is one of the first to be observed and 
published with four JWST observing modes: NIRSpec (Bóker 
et al. 2023), NIRCam (Rieke et al. 2023), MIRI (Wright et. 
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Table 1. Source Information for CEERS. 1019 


R.A. 215.0353914 [deg] (1) 
Dec. 52.8906618 [deg] Q) 
mE160W 25.2 + 0.03 [AB mag] (3) 
mE356W 24.8 + 0.01 [AB mag] (4) 
Zphot 8.844012 (HST+IRAC) (5) 
Zphot Bg (+JWST) (6) 
Zspec(Lyo) 8.683+0:001 (Keck) (7) 
Zspec(Lya) 8.6854 + 0.0045 (NIRSpec) (8) 
Zspec([O m] 8.6788 + 0.0002 (NIRSpec) (9) 
log(Mx) 9.5 + 0.3 [Ms] (10) 
log(sSFR) -7.9 + 03 [Mo yr-!] (11) 
log(Mgy) 6.95 + 0.37 [Mo] (12) 
Te({O m]) 18630.76 + 3.68 [K] (13) 
12+log(O/H) 7.66 + 0.51 (14) 
Z 0.095+0:21 [Zo] (15) 
Ne 1.9+0.2 x 10? [em73] (16) 
UV Slope £ -1.76:0- 2 (17) 
Ay 0.4 € 0.2 [mag] (18) 


NorE— (1) Right ascension. (2) Declination. (3) Observed AB 
magnitude in the F160W filter from Finkelstein et al. (2022a). (4) 
Observed AB magnitude in the F356W filter (82.2). (5) Photometric 
redshift measured with HST (prior to JWST photometry) Finkelstein 
et al. (2022a). (6) Photometric redshift measured including JWST 
photometry from CEERS (§5). (7) Spectroscopic redshift measured 
from Keck/MOSFIRE spectroscopy via Lya emission line detection 
from Zitrin et al. (2015). (8d) Spectroscopic redshift measured from 
JWST/NIRSpec spectroscopy via Lya. (9) Spectroscopic redshift 
measured from JWST/NIRSpec spectroscopy via [O m] emission 
line detection (§4.1). (10) Stellar Mass from Prospector SED fit 
(85). (11) Specific star formation rate ($5). (12) Black Hole Mass 
($6.1). (13) Electron Temperature (86.2). (14) Te-based Metallicity 
($6.2). (15) Metalliticy ($6). (16) Electron density ($6.2). (17) UV 
Spectral Slope ($5.1). (18) Dust attenuation (85 & 6.4). 


al 2023, submitted), and NIRCam/WFSS. We describe these 
observations below and provide a summary of information 
about this source in Table 1. 

Additional imaging data for this source was obtained in 
the rest-frame infrared from Spitzer/MIPS at 24um and Her- 
schel/PACS at 100 um, as well as SCUBA-2 at 850um and 
the JVLA at 3GHz which is shown in Appendix A. X-ray 
imaging from the Chandra Space Observatory for this source 
is also discussed in Appendix B. 


2.1. NIRSpec Observations 


The source presented in this work is included in the 
JWST/NIRSpec (Jakobsen et al. 2022) MOS configurations 
taken with the Micro Shutter Array (MSA; Ferruit et al. 2022) 
during the CEERS epoch 2 observations (December 2022). 


These NIRSpec observations are split into 6 different MSA 
pointings, each of them observed with the G140M/F100LP, 
G235M/F170LP, and G395M/F290LP medium resolution 
(R « 1000; here denoted by “M”) gratings plus the prism 
(R ~ 30-300), fully covering the ~ 1-5 um wavelength 
range. The MSA was configured to use 3-shutter slitlets, en- 
abling a 3-point nodding pattern, shifting the pointing by a 
shutter length plus the size of the bar between the shutters in 
each direction for background subtraction. The total exposure 
time per disperser is 3107s distributed as three integrations 
(one per nod) of 14 groups each in the NRSIRS2 readout 
mode. 

JWST/NIRSpec 2D+1D spectra for this source in each of the 
M gratings (G140M, G235M, G395M) are shown in Figure 1 
with the location of the source marked by two horizontal red 
lines in each 2D spectrum (this source was not observed with 
the prism). Details of the extraction method from 2D to 1D 
are detailed in $2.1.1 below. 


2.1.1. NIRSpec Data Reduction 


The details of the CEERS NIRSpec data processing will be 
presented by Arrabal Haro et al., (in prep.). We summarize 
the main steps of the reduction here. The NIRSpec data is 
processed with the STScI Calibration Pipeline! version 1.8.5 
and the Calibration Reference Data System (CRDS) mapping 
1027. We use the calwebb_detector1 pipeline module to 
subtract the bias and the dark current, correct the 1/f noise 
and generate count-rate maps (CRMs) from the uncalibrated 
images. At this stage, the parameters of the jump step are 
modified for an improved correction of the “snowball” events? 
associated with high-energy cosmic rays. 

The resulting CRMs are then processed with the 
calwebb. spec2 pipeline module, which creates two- 
dimensional (2D) cutouts of the slitlets, performs the back- 
ground subtraction making use of the 3-nod pattern, corrects 
the flat-fields, implements the wavelength and photometric 
calibrations and re-samples the 2D spectra to correct the dis- 
tortion of the spectral trace. The pathloss step accounting 
for the slit loss correction is turned off at this stage of the re- 
duction process (see Figure 2). Instead, we introduce slit loss 
corrections based on the morphology of the sources in the 
NIRCam bands and the location of the slitlet hosting them. 

The one-dimensional (1D) spectra of the sources were ob- 
tained via an optimal extraction (Horne 1986) with a spatial 
weight profile from the trace of the source such that the pixels 
near the peak of the trace are maximally weighted. To create 
the extraction profile, the 2D signal-to-noise (SNR) spectrum 
was collapsed in the spectral direction for each grating in- 


! https://jwst-pipeline.readthedocs.io/en/latest/index.html 


? https://jwst-docs.stsci.edu/data-artifacts- and- features/ 
snowballs-and-shower-artifacts 
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Figure 1. 2D and 1D spectra of CEERS. 1019 from three /WST/NIRSpec M gratings, G140M (top), G235M (middle), and G395M (bottom). 
The horizontal red dashed line identifies the central location of the source in the 2D spectrum and is the extraction center for our 1D spectra. 
Description of the CEERS NIRSpec observations for this source are in $2.1, and the data reduction process is described in §2.1.1. 


dependently, taking the median value at each spectral pixel 
and fitting a Gaussian to the positive trace. This source is ef- 
fectively unresolved at all wavelengths (intrinsic size FWHM 
is less than the JWST PSF FWHM), so the central trace is 
unaffected by upper and lower negative traces. A significant 
trace was measured in both the G140M and G395M gratings 
but was not as apparent in the G235M such that it could be 
fit with a similar Gaussian profile. We thus used the profile 
from the G140M grating for the extraction from 2D to 1D in 
the G235M grating. 

We test the accuracy of the error spectrum computed by the 
pipeline by comparing the normalized median absolute de- 
viation of the science spectrum (masking out emission lines, 
and removing a smoothed continuum) to the median of the 
error spectrum in each grating individually. The real data 
show fluctuations ~1.5—2x larger than the typical error value. 
Thus, we measure and scale the error spectrum up by this 
scale factor in each grating. 


We evaluate the scaling factor for the slit-loss correction 
compared to the NIRCam photometry for this source, again 
using the reduction without the path loss correction applied. 
In all NIRCam filters, we repeatedly calculate the fluxes en- 
closed within the 072 x 0"46 rectangle aperture using the 
source and MSA shutter positions and estimate the scaling 
factors to match them with the total flux measurements. The 
scaling factors range from approximately 2.0 to 2.5 among 
the NIRCam filters. We correct the spectrum with the scaling 
factor of the filter whose central wavelength is closest to that 
of the observed wavelength. This self-consistently applies an 
aperture correction accounting for the variable PSF across the 
observed wavelength range. 

Once the errors are corrected, and the spectra are scaled to 
the NIRCam photometry, the three M gratings are combined 
into a single spectrum, re-sampling to a common wavelength 
array at the overlapping wavelengths and adopting the mean 
flux at each pixel, weighted by the flux errors. This combined 
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Figure 2. Comparison of the JWST/NIRSpec combined M grating spectrum when using the native Pipeline pathloss correction (black) vs. 
without this step of the pipeline but rather scaling the spectra to the measured JWST/NIRCam photometry (blue) as described in Fujimoto 
et al. (2023) and in §2.1.1. Our scaled spectrum (purple) highlights how default pipeline pathloss correction under-predicts the total slit-loss 
correction for this source and that corrections to the NIRSpec spectrum are required for flux calibration of resolved sources. 


spectrum is then used for the remainder of the analysis in this 
paper. 


2.2. NIRCam Imaging 


The galaxy discussed in this work was observed in the 
CEERS JWST/NIRCam (Rieke et al. 2003, 2005; Beichman 
et al. 2012; Rieke et al. 2023) imaging taken in December 
2022 (in CEERS NIRCam Field 8); these imaging data, in- 
cluding the detailed reduction steps, are described in Bagley 
et al. (2022). Photometry of this source was measured us- 
ing Source Extractor (Bertin & Arnouts 1996) in dual-image 
mode. The photometry procedure is broadly similar to that 
described in full by Finkelstein et al. (2022b), though with a 
few differences we detail here. First, to improve color accu- 
racy, especially in the HST/WFC3 bands, we do a two-step 
point-spread function (PSF) correction. For any filters with a 
PSF full-width half max (FWHM) smaller than the NIRCam 
F277W (which includes ACS F606W, F814W, and NIRCam 
F115W, F150W and F200W), we convolve the images with a 
kernel designed to match these images’ PSFs to that in F277W. 
For images with larger PSFs (WFC3 F105W, F125W, F140W 
and F160W, and NIRCam F356W, F410M and F444W), we 
derive a correction factor. This is done by convolving the 
F277W image with a kernel designed to match the PSF to 
the larger PSF image, with a correction factor derived as the 
ratio of the flux in the broadened image to that in the original 
F277W image. In this way, we correct for light not captured 
in our default aperture in the bands with larger PSFs, with- 
out needing to smooth all bands to that with the largest PSF 
(which would be WFC3 F160W). 

We tested this process using source injection simulations 
to confirm that accurate colors were recovered. Finally, fol- 
lowing Finkelstein et al. (2022b), we derived accurate total 
fluxes first by correcting the fluxes in all bands by an aperture 


correction derived in F277W (as the ratio between the flux 
in our small Kron aperture used to measure colors, and in 
the default larger "FLUX. AUTO" Kron aperture), and then 
correcting this by a factor of 8—10% to account for missing 
light from the wings of the PSF, with this correction factor 
derived via the source injection simulations (see Finkelstein 
et al. 2022b for more details). 

This final photometry catalog includes measurements over 
the full CEERS NIRCam wavelength range in the F115W, 
F150W, F200W, F277W, F356W, F410M, and F444W filters, 
which have exposures of ~ 3000 seconds per filter (~ 6000 s 
for F115W); as well as in the existing HST/CANDELS ACS 
and WFC3 F606W, F814W, F105W, F125W, F140W, and 
F160W bands. We include the F098M (non-detection) data 
for this source from Finkelstein et al. (2022a). 

The measured photometric redshift for this galaxy before 
JWST was Zpnhor (HST) = 8.84*07? (Finkelstein et al. 20222). 
With the addition of the JWS T/NI IRCam i imaging (and remov- 
ing the blended Spitzer/IRAC data), we measure a photomet- 
ric redshift with EAZY (using the same process as Finkel- 
stein et al. 2022b) nearly equal to the spectroscopic redshift 


(Zspec = 8.679) of Zphor (HST + JWST) = 8.68:002. 


2.3. MIRI Imaging 


This source was also observed with JWST/MIRI in CEERS 
Epoch 1 with filters F560W and F770W. The photometry 
is obtained from the CEERS team DRO.5 data release for 
the MIRI 3 and MIRI 6 fields (G. Yang et al. in prep)>. 
The MIRI data reduction process is described by Yang et al. 
(2021), and the photometric measurements for this source 
are previously reported by Papovich et al. (2022). We 


3 https://ceers.github.io/releases.html 
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Instrument Filter Flux (nJy) 

HSTIACS F606W 8.95 + 7.67 
F814W 0.96 + 8.07 

HSTIWFC3 F098M -11.00 + 18.40 


F125W 200.78 + 11.30 
F140W 272.05 + 19.18 
F160W 306.22 x 9.69 
JWST/NIRCam F115W 117.07 + 3.81 

F150W 334.15 € 7.11 

F200W 313.72 + 5.63 

F277W 371.35 € 4.51 

F356W 449.55 € 3.51 

F410M 450.11 + 8.81 

F444W 1039.16 + 7.65 
JWST/MIRI F560W 426.80 + 21.60 
F770W 404.00 + 16.90 

Table 2. Photometric measurements for CEERS_ 1019, with fluxes 
in nJy. The HST photometry is updated from Finkelstein et al. 
(2022a) using the NIRCam-selected apertures as described in Finkel- 
stein et al. (2022b). JWST/NIRCam photometry from CEERS Epoch 
2 imaging measured in a similar way as Finkelstein et al. (2022b), 
described in $2.2. JWST/MIRI photometry from Papovich et al. 
(2022) using CEERS Epoch 1 imaging, described in §2.3. 


present the JWST/NIRCam and MIRI photometry for this 
source in Table 2. With the addition of the MIRI pho- 
tometry, the photometric redshift is slightly modified to 


Zphot (HST + JWST) = Bu. 


2.4. NIRCam Grism Observations 


This source was also observed with the JWST/NIRCam 
Wide-Field Slitless-Spectrograph (WFSS) in December 2022 
as part of the CEERS Epoch 2 data. These data were obtained 
using the F356W filter with both the orthogonal column “C” 
and row “R” grisms. Grism spectra for this source are shown 
in Figure 3; the data reduction methodology will be presented 
in N. Pirzkal et al. (in prep). We observe several emission 
lines detected in the combined grism spectrum. These are 
consistent with [O II], [Ne III], He and Hó if a wavelength 
shift (of ~60A) is applied from the NIRCam to NIRSpec 
spectra. This is likely due to as-yet improving wavelength 
calibrations for both instruments, though given the agree- 
ment of our observed Lya wavelength with ground-based 
measurements, the offset is likely dominated by the NIRCam 
WESS spectra. While upcoming improved calibrations will 
provide the necessary adjustments to compare better the re- 
sulting line fluxes and detections from JWST, these detected 
lines show the utility of this NIRCam WFSS mode for early- 
galaxy spectroscopy. For the remainder of this paper, we use 
all measurements from the higher signal-to-noise NIRSpec 
spectra. 


3. METHODS: EMISSION LINE SEARCH 


To search for emission line features in our 1D spectrum, 
we utilize an automated line-finding code first published in 
Larson et al. (2018) and outlined here. This code uses a 
Monte Carlo Markov Chain (MCMC) routine to fit a model 
which consists of a Gaussian line plus a continuum constant 
to a given wavelength range, with four free parameters: the 
continuum level, line central wavelength, line FWHM, and 
integrated line flux. To run the MCMC, we use an IDL 
implementation of the affine-invariant sampler (Goodman & 
Weare 2010) to sample the posterior similar to that used in 
Finkelstein et al. (2019), which is similar to the Python emcee 
package (Foreman-Mackey et al. 2013). 

As a first pass-through for emission line features, we search 
the entire spectrum with our automated code. At each wave- 
length pixel, we do an initial SNR check (flux/error). If this is 
greater than unity, we run our fitting routine (this "pre-check" 
is not required but helps with efficiency as no detectable emis- 
sion feature would have SNR < | at the line center). We use a 
fitting range of 100x the pixel scale in the G140M grating on 
either side of the center pixel, which equates to 630 À across 
the full spectrum. This ensures we fit the same wavelength 
range for every line, regardless of pixel scale. We fit a single 
Gaussian feature to the spectrum: 


_ 442 
1 an (1) 


f(a) = fe+ foexp (- 2 G2 


We also impose non-informative priors on each free pa- 
rameter designed to increase the efficiency of the line de- 
tections while removing the chance for the MCMC chain 
to exit a range of realistic parameters. For the continuum 
constant (fe), we let it vary between —3 to 10x the average 
flux across the fitting range. These values are larger than 
the 10 noise level and broadly encompass the typical contin- 
uum values for our source. We restrict the peak wavelength 
(Ao) to be the wavelength at that pixel + one pixel, such 
that we fit a Gaussian within each pixel. The actual value 
in A varies as the pixel scale (AA) for NIRSpec is wave- 
length dependent with AAgi4om = 6.4 A, AAqoasM = 10.7 
A, and Adg3osm = 17.9 A. We limit the full-width half-max 
(FWHM = 2,/21n(2)o = 2.3550) to +30 km s^! of the ob- 
served FWHM of the [O u1]500g line, the brightest feature in 
the spectrum (see Figure 4 and Table 3). Tying the FWHM 
to [O iii]soos in velocity units rather than pixels is critical as 
the wavelength dispersion changes across the observed wave- 
length range. The line flux prior requires the line flux (fo) 
to be greater than the average of the flux over the plotting 
range after a > 1c clipping (estimate of the continuum at that 
location) and less than 100x the average flux, well outside the 
expected range of line fluxes. 

To define the starting point for our MCMC routine, we used 
the IDL routine mpfit, a Levenberg-Marquardt least-squares 
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Figure 3. The spectrum of CEERS_1019 obtained with JWST/NIRCam’s wide-field slitless spectroscopy mode. This consists of spectra with 
both the column (C; blue) and row (R; purple) grisms, both taken with the F356W blocking filter, with a combined spectrum shown in black. 
While these data have a higher noise level than the NIRSpec spectra, we observe detections of the same [O 11]3727+3729, [Ne 111]3g79 and bluer 
Balmer lines that we see in NIRSpec (Figure 5). While we do not use the grism measurements in this paper, these data highlight the utility of 
this mode for obtaining slitless measurements of modestly faint emission lines out to high redshifts. 


fitting routine that fits a Gaussian with the above parameters 
(Markwardt 2009). We then run our MCMC fitting code with 
10,000 iterations and 100 walkers on each pixel and determine 
the best line fit results as outlined in the following sections. 
We use the median of the last 10,000 steps of our MCMC 
chain for our fit parameters. To measure the MCMC error on 
our parameters, we use the robust sigma calculation: using 
the median absolute deviation as the initial estimate, then 
weighting points using Tukey's Biweight (equation 9 from 
Beers et al. 1990), assuming a Gaussian posterior distribution. 

We place several constraints on the resulting Gaussian fits 
when determining a successful line detection, including the 
wavelength and comparison to neighboring pixels. First, we 
mask out the edges of the wavelength range for the three M 
Gratings (within 0.97 — 5.1m), where the grating transmis- 
sion curve falls below 70%. To determine the correct peak 
pixel for the line, we select the fit at the pixel which has the 
highest peak flux (flux at the peak of the Gaussian fit) within 
the FWHM of the line. We calculate an integrated SNR from 
the MCMC Gaussian fit as the median line flux divided by 
the robust sigma of the line flux as described above. We also 
measure a peak SNR value, calculated by taking the maxi- 
mum flux of the Gaussian fit and calculating the ratio of this to 
the noise per pixel within our fitting range (similar to Larson 
et al. 2022). For our initial search through the spectrum for 
emission lines, we only consider lines found at > 5c in both 
integrated and peak SNR. This initial pass is run for mask- 
ing out significant features in the spectrum to fit and remove 
the detected continuum, as described in $3.1, and identifying 
potentially unexpected emission line detections in sources at 
this redshift. In $3.2, we discuss how these initial SNR mea- 
surements overestimated the significance of our emission line 
detections and how we determine final line flux errors and 
SNR values for our detected emission lines. 


3.1. Measuring Continuum in NIRSpec 


To fit the continuum to the NIRSpec spectra, we first mask 
out all the emission lines (peak of the line + FHWM) and 
interpolate over those pixels with an average of the 3 not- 
masked pixels on either side of the line region. We then 
smooth the entire spectrum using a boxcar filter with a width 
of 60 pixels. As each M grating has a different pixel scale, 
we then smooth each section of the spectrum with a larger 
boxcar filter at the blue end compared to the red end. We 
use 120 pixels where 2 < 1.75um (the G140M grating), 60 
pixels from 1.75 < 4 < 2.9um (the G235M grating), and 30 
pixels where A > 2.9um (the G395 grating). To remove any 
discontinuous jumps between the gratings, we smooth the 
whole spectrum again with the 60-pixel boxcar filter. This 
measured continuum is plotted in green over the combined 
NIRSpec spectrum in Figure 4 with notable detection in the 
G140M grating but non-significant detection in the redder 
gratings. 


32. Determining the Significance of the Emission Lines 


Upon inspection of the initially identified lines described 
above, it was found that this methodology tended to over- 
estimate the signal-to-noise ratio (SNR), especially for faint 
lines. We thus empirically derive line flux errors and appro- 
priate SNR measurements for each line that passes the above 
selection cuts in the following way. Using the line-masked 
spectrum as above, we find the 20 closest "blank" pixels, free 
of any other nearby detected emission lines, on either side 
of the emission line (line center + FWHM such that we are 
not overlapping the emission line at all). At each of these 
blank pixels, we insert a fake emission line with the same 
parameters as our detected emission line (F(A, FHWM, and 
Jc) and a line-center (Ao) at the wavelength of that pixel. We 
then run the same line-fitting routine at each pixel and record 
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Figure 4. Combined JWST/NIRSpec spectrum from G140M+G235M+G395M of CEERS. 1019, plotted in F} [1071? erg s7! cm7? A7!] vs 
Observed Wavelength [um]. Plotted in green is the fit to the continuum in the spectrum (see $3.1) plus the detected emission lines as discussed 
in $4. We require a SNR >2.3 (with our simulation-estimated noise as described in $3.2) for a line to be considered detected. 


the recovered line flux from the MCMC (7F(A)54) for each of 
our 40 simulated lines. Our reported line flux error (AF(4)) 
is then the median absolute deviation of these recovered line 
fluxes, as reported in column 4 of Table 3. The SNR of 
our emission lines is taken as the measured line flux from 
our line fit divided by this error, SNR = F(A)/Af(A) and is 
reported in column 3 of Table 3. To derive a consistent thresh- 
old in SNR where we consider a line detected, we examine 
the results of the automated line-finding routine, searching 
for lines with rest-frame wavelengths which do not match a 
known line (using a large line list). The point below which we 
detect no “unknown" (spurious) emission lines is SNR = 2.3, 
which we use as our detection threshold for emission lines in 
CEERS 1019 unless otherwise noted below. 


3.3. Rest-frame Equivalent Width Measurements 


To measure rest-frame equivalent widths for our emission 
lines (EW est = za. we require a measurement of the 
continuum emission in the region near a given emission line. 
While, in principle, this could come from the spectrum it- 
self, the continuum is only detected at a high significance in 
the G140M (blue) grating. Thus, to enable a uniform pro- 
cedure for all emission lines, we use the best-fit Prospector 
model (see $5 and Figure 7) for our rest-frame equivalent 
width measurements. We note that this is based on the same 
NIRCam photometry used to derive our slit-loss corrections, 
and we also find that the agreement with the spectroscopic 
continuum in G140M is excellent. To determine the contin- 
uum value of the Prospector model, we first mask out any 
emission lines and smooth the spectrum following the same 
process described above. Our EW continuum value is the 
median of the Prospector continuum over the width of the 
line. The error is taken as the median absolute deviation of 
the Prospector model over this same range. 


3.4. Non-Gaussian Emission Line Fits 


In our automated line search above, we fit a single Gaussian 
to each emission line, which is not physical for many of the 
lines in our spectrum. We thus implement three additional line 
profiles: a Dual-Component Gaussian, a Doublet Gaussian, 
and an Asymmetric Gaussian fit as described below. The fits 
are done using the same MCMC routine as above, and any 
parameters not mentioned below as having been changed or 
added are the same as those in the single Gaussian fits. 


3.4.1. Asymmetric Gaussian Fit 


Given the detection of Lya in our NIRSpec spectra, we 
want to fit an appropriate Asymmetric Gaussian profile to the 
data to account for IGM absorption on the blue side of the 
line. For this fit, we add an extra FWHM parameter such that 
we now have a ‘blue’ and a ‘red’ side FWHM of the line. 
The blue FWHM can vary between 0, as the blue side of the 
line may be fully attenuated, and FWHM) m] + 30 km s7!, 
which is the same maximum as the single-Gaussian fit. We 
start the MCMC chain for this parameter at the FWHM) ui]: 
The red FWHM is less restricted as it can vary between the 
FWHM/9 ur] and 2000 km s~!, and the MCMC chain is 
started at FWHMr[o ui] + 30 km s~!, the maximum value for 
the blue FWHM. Following the same method as Jung et al. 
(2020) , the equation for the Asymmetric Gaussian is: 


exp ( - m: for A < Ao, 
f(A) asym = fe + fo X (ata)? (2) 
exp(- A ) for A > 4o, 
red 


This returns the continuum flux (fe), peak line flux value 
(fo), a peak wavelength (40), and the blue- and red-side line 
widths (ciue and Crea). To measure our line flux, we integrate 
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Table 3. Measured Line Fluxes in CEERS_1019 


Line Aobs SNR Line Flux FWHM Flux Flux FWHM Rest EW Av 
Arest Central/Blue Total Single/Blue Narrow/Blue Broad/Red Broad/Red 
[Å] [107!8 cgs] [km s71] [10718 cgs] [10718 cgs] [km s71] [Å] [km s71] 
(1) (2) (3) (4) (5) (6) (7) (8) (9) (10) 
Lya 1215 11774.30 + 5.49 3.70 5.44 + 0.03 207.88 + 111.50 - - 1360.66 + 479.93 10.73 € 0.56 206.20+139.97 
N V123841242 11969.40 + 3.39 1.98 1.742 0.88 358.93 + 19.66 0.03 €0.03 1.70 € 0.72 - 3.90 + 2.30 -522.45+85.20 
Niv]i48341486 14357.60 + 8.77 7.59 3.36 € 0.44. 363.99 + 20.77 1.23 +0.85 2.13 + 0.69 - 7.96 + 1.47  16.15x183.26 
C Il] 59g 18493.30 + 9.22 1.63 1.702 1.05 642.90 + 156.84 - - - 6.04 € 3.88 309.514149.63 
Mg n2796+2803 27077.90 + 6.67 1.18 1.62 =+ 1.37 358.76 £20.82 0.67+1.01 0.96 + 1.20 - 11.56€ 9.90 144.74€74.18 
[O 11]372743729 36074.80 + 4.09 11.86 2.94 +0.25 357.55 + 20.83 1.80+0.35 1.15 € 0.35 - 34.79 + 5.80 11.43434.71 
[Ne i1]3869 37456.30 + 2.29 12.84 3.48 +0.27 358.44 + 17.43 - - - 44.27 x 7.07 8.51419.64 
He 13gg9+H7 38909 37656.80 + 9.73 2.28 0.66+0.29 357.13 + 18.12 - - - 8.43 + 3.87 . 71.95x77.78 
[Ne 11]396g+He397; 38407.90 + 3.75 5.48 2.04€ 0.37 347.65 21.00 0.48 +0.43 1.56 € 0.40 - 26.53 + 6.07 -23.34+30.11 
Hó4102 39717.90 + 6.07 3.06 1.14 € 0.37 358.56 + 18.29 - - - 15.64 + 5.50  53.32+46.36 
Hy4341 42018.20 + 5.25 6.59 2.60+0.39 360.56 + 18.78 - - - 39.36 + 848 -27.63438.12 
[O 11] 4364 42244.90 + 9.91 3.03 1.16+0.38 361.11 + 17.79 - - - 17.78 +6.45 19.06-70.68 
HB4862 47059.80 + 4.91 10.21 7.32 +0.72 354.80+19.61 3.79 €0.69 3.53 € 1.51 1196.26 + 349.07 128.94 + 18.96 -31.55+32.06 
[O m1] 4960 48009.60 + 1.87 26.34 13.08 + 0.50 358.52 + 17.37 - - - 235.35 + 27.51 1.00+13.64 
[O ir1]500g 48473.50 + 1.14 47.16 39.58 + 0.84 358.66 + 12.23 - - - 725.59 + 82.21 -0.59 + 124.41 


Note— Identified emission line values in CEERS_1019, see Figure 5 for individual plots. (1) Line name and rest-frame vacuum wavelength in 
[A]. (2) Observed wavelength of the line in NIRSpec in [A]. (3) Signal to Noise (SNR) of the line (§3.2). (4) Total integrated line flux in Fy 
units of [10718 erg s! cm7? A-!] ($3). (5) FHWM of the line fit in [km s^! ] for single-Gaussian fits, the blue line in the Doublet (83.4.2) and 
Asymmetric (§3.4.1) fits, or the narrow component in the Dual-Component (§3.4.3) fits. (6) Flux of the line fit in [107!8 erg sem? À-] 
for single-Gaussian fits, the blue line in the Doublet ($3.4.2) and Asymmetric (§3.4.1) fits, or the narrow component in the Dual-Component 
(83.4.3) fits. (7) Flux of the red line in the Doublet (83.4.2) fit or the broad component in the Dual-Component (83.4.3) fit in [10718 erg s7! 
cm~? wi (8) FHWM of the red line in the Asymmetric (83.4.1) fit, or the broad component in the Dual-Component (83.4.3) fit in [km sl]. 
(9) Rest-frame Equivalent Width (Rest EW) of the line in [A] (83.3). (10) Velocity offset of the emission line from the systemic redshift of 
CEERS 1019 (zt m] = 8.679, 84.1) in [km s^! ]. 
*While the C mı] line is a Doublet with peaks at 1906 and 1908 À in the rest frame, we do not resolve this doublet and show the best fit, which 
uses a single-Gaussian profile with a broad FWHWM ($84). 


this emission line profile and report the integrated line flux in 
column 4 of Table 3 


3.4.2. Doublet Gaussian Fits 


Several of the UV and Optical emission lines are not a single 


line flux for the doublet (f(A)) as reported in column 4 in 
Table 3, as well as individual line fluxes for each component, 
as reported in columns 6 and 7 of that same Table. 


3.4.3. Dual-Component Gaussian Fit 


emission line feature but rather a doublet. To accurately fit 
these lines, we implement a doublet Gaussian line profile to 
our MCMC routine: 


= 2 
f (A)Doublet = fc + fi exp ( E ; A -9 ) 
1 (4 — (41 + 6A))* e 
*Aee(-57— — —) 


For this fit, we employ the continuum flux (f) plus two 
Gaussians with the same FWHM (o^, individual line fluxes 
for each line (fı and f» respectively), and with the blue peak 
(41) fixed to +6 pixels from the expected location, and the red 
peak fixed to the redshifted separation between the two dou- 
blet lines (6A = Arest(1 + z)). This returns a total combined 


Upon inspection of the initial emission lines, the Hf line 
had a noticeable second broad component feature, as dis- 
cussed in $4.2 below. In order to accurately measure this, we 
employ a dual-component Gaussian profile: 


1 (A = 4g? 
f G)pua = fe + fnac EXP ( — ula) 
O narrow 
2 (4) 
1 (A—- Ag) 
* Soroad exp ( m 29 
Ü broad 


This fit utilizes the continuum flux ( f.) plus two Gaussians 
with the same line center (Ag), individual line fluxes for each 
line component (fnar and foroaa respectively), and separate 
FWHM for each (Onar and Oproad). The narrow FWHM has the 
same constraints as the single Gaussian fit (+ FHWM/9 ui 
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but the broad component is not allowed to be smaller than 
the narrow component (>FHWM/o9 ui] * 30 km s^!) and can 
extend to the unrestricted FWHM (< 2000 km s^! ). We start 
the MCMC chain for the narrow FWHM at FHWM/9 ui]: 
and the FWHM of the broad component at 3x this value. 
This returns a total combined line flux for the line (F(A)) as 
reported in column 4 in Table 3, as well as individual line 
fluxes for each component, as reported in columns 6 and 7 of 
that same Table. 


4. EMISSION LINE MEASUREMENTS 


We find abundant emission lines in CEERS 1019 and 
present our measured line fluxes in Table 3. A plot indicating 
the detected emission lines is shown in Figure 4, and fits to 
each line are shown in Figure 5. In the following sections, we 
briefly discuss the fit to different lines and provide a detailed 
description (where necessary) of the detected emission lines. 
We note that H£ has a significant broad component, which 
we discuss in detail in $4.2. 


4.1. Redshift Confirmation 


We determine the spectroscopic redshift using the fit 
to the [O ]soos emission line, the brightest in our spec- 
trum. This line is detected with a peak wavelength of 
A = 48473.50 + 1.14 A, yielding a spectroscopic redshift 
for this source of zio yy] = 8.6788 + 0.0002. When compar- 
ing this to the redshift we measure from the line centers of 
the other strong [O m]4960 line we get a consistent measure- 
ment within Az/z < 107^. We use this value as the systemic 
redshift for CEERS_1019. 


4.2. Broad Hf Feature 


The observed profile of the Hf emission line shows a strong 
narrow component, similar to [O 11], but also a weaker broad 
component. To explore the significance of this broad compo- 
nent, we perform a dual Gaussian fit to this line as described 
in $3.4.3. Uncertainties on both components were derived 
via similar simulations as described in $3.2. From this fit, we 
find a distinct broad component to the Hf emission line, with 
a FWHM = 1196 + 349 km s^!. We show this fit (green) in 
Figure 6. The flux in this broad component is comparable 
to the flux in the narrow component, with a broad-to-narrow 
flux ratio of 0.93. 

We consider multiple measures of the significance of this 
broad feature. The first is the signal-to-noise of the broad line 
flux from this fit, which is 2.5. However, this does not con- 
sider how good the fit is without the broad line component. 
We next try a narrow-line only fit with the FWHM allowed to 
be free (i.e., not tied to [O iri]) shown as the blue line in Figure 
6. We find a single-line FWHM of 500 + 52 km s^, signif- 
icantly broader than the fiducial narrow-line-only FWHM of 
367 + 17 km s^! (purple line in Figure 6), making it clear 


that the broad component affects even a single-line fit. Fi- 
nally, we calculate the Bayesian Information Criterion (BIC) 
between the dual-component and narrow-only line fits. The 
BIC is a method of comparing the goodness-of-fit between 
two models, accounting for differing degrees of freedom. It 
is defined as BIC = x? + k In(N), where k is the number of 
free parameters, and N is the number of data points (number 
of spectral pixels in this case; Liddle 2004). Our narrow-only 
fit with four free parameters (fe, fo, 4o, and o - see $3) has a 
BIC of 21.4 when the FWHM is tied to [O m] and 21.8 when 
the FWHM is free. The narrow+broad fit with six free pa- 
rameters (fe, Ao, fnar» foroads C nar; ANd Ohroad - See $3.4.3) has 
a BIC of 18.7. The ABIC = 2.7-3.1 gives "positive" evidence 
(using the scale of Jeffreys 1961) that a broad-line component 
is necessary to fit the Hf line successfully. 

This significant broad component of Hf thus indicates the 
presence of high-velocity gas, which we interpret as emitting 
from the broad-line region of an AGN. Large-scale outflows 
can similarly produce broad velocity widths but would be ap- 
parent in all emission lines (e.g. Amorín et al. 2012; Hogarth 
et al. 2020). We do not observe any broad features in lines 
with much higher signal-to-noise like [O 11]A5008. Instead, 
[O m] and other forbidden lines have widths that are incon- 
sistent with the outflow scenario and are instead consistent 
with the narrow Hf component. Broad-line AGN typically 
exhibit broad components for permitted lines (like H£) and 
narrow widths of forbidden lines (like [O mJ) (e.g. Schmidt 
et al. 2016; Vanden Berk et al. 2001), as observed in our 
spectrum. We discuss further evidence for the AGN nature 
of this source in $ 4.4 below. 


4.3. Lya 


The Lya line was first detected by Zitrin et al. (2015) at 
ZLya = 8.683 with Keck/MOSFIRE. We detect an emission 
line at this same wavelength, which we fit with an asymmetric 
Gaussian profile as described in §3.4.1, and measure a line 
flux ÁLya = 5.44 + 0.03 x 10° erg s^! cm ?. Using our 
measured central wavelength and the vaccum wavelength for 
Lya of 1215.67 A, we measure a Lya-based redshift of 8.6854 
+ 0.0045, 260 + 143 km s^! higher than that of [O ni]. This is 
lower, although consistent within ~1o of the measured a Lya 
velocity offset (relative to N v) of +362 km s7! from Mainali 
et al. (2018). 

We find a clear asymmetric line profile with an extended 
red-side tail (FWHM,eq> 1000 km s~!). Although it becomes 
more difficult to detect Lya at this high redshift due to the in- 
creasingly neutral IGM, one can expect the escape of Lya that 
is scattered and redshifted enough to avoid resonant scattering 
(Dijkstra et al. 2004), as seen in the extended red tail in our 
Lya spectra. Additionally, a sharp cutoff at the blue edge of 
the line profile may indicate a significant contribution of reso- 
nant scattering, which can be caused by a proximate optically 
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Figure 5. Fits to the emission lines identified in the JWST/NIRSpec spectrum of CEERS_1019. Each panel shows an individual emission line 


fit, with the type of line profile as described in §3. The panels are plotted in F} [10~!9 erg s7 


1 


cm ? À- vs Observed Wavelength [um], are 


presented in order of increasing wavelength, and are scaled vertically to show the extent of the highlighted emission line. Emission line values 
are tabulated in Table 3 and discussion of specific lines is in $4. 
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Dual Component 
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Figure 6. Left: | Dual-component Gaussian fit to the Hf emission line (green) yielding a FWHMpyoaq = 1196.3 + 349.1 km s^! and 
FWHMpar = 354.8x 19.6 km s-l. Alternative single-Gaussian line profile fits are shown where the FWHM is restricted to FWHM m] = 30 


km s^! (blue) or left unrestricted (purple). In both cases, the fit to the H8 emission line is worse than the dual-component fits as indicated by 
the higher Bayesian Information Criterion (BIC) measurement (§4.2). Right: Scaled down fits of the same dual-component fit from Hf to the 
peak of the Hó (top) and Hy (bottom) emission lines. In both cases, the broad component is not distinguishable from the noise, indicating that 
a broad component could be present in these lines below our noise threshold. 


thick medium (Mason & Gronke 2020). Particularly, reion- 
ization simulations predict the location of the sharp blue-side 
edge could be redshifted due to the infall motion of neutral 
gas around a galaxy, as seen in our Lya spectrum. The cutoff 
location is found at ~200 km s~! in our spectra, which is com- 
parable to the predicted velocity of the in-falling gas around 
a Muy = —22 galaxy at this redshift (Park et al. 2021). 

Our Lya line flux is ~7 times fainter than the line flux 
measured from MOSFIRE. This could be due, in part, to the 
MOSFIRE slits being wider than the NIRSpec micro-shutters 
(0.7" vs. 0.2"). Extended Lya emission in this source would 
not be corrected for by our NIRSpec slit-loss correction. The 
MOSFIRE observations of this source were taken ~ 8 — 9 
years prior to the JWST/NIRSpec data, ~ 1 year in the rest- 
frame, and it is possible that the variability common in AGN 
might contribute to the difference in the measured line fluxes. 
AGN variability is typically ~20% (e.g., MacLeod et al. 2012) 
and is unlikely to fully explain this large difference. 

The Nvj243 line was detected from this source with 
Keck/MOSFIRE by Mainali et al. (2018) at 12019.5 A, with 
the other line in the doublet, N v1238, obscured by a skyline 
at 11981A. In our JWST/NIRSpec data, we do not detect a 
line at either of these wavelengths, but there is a 1.80 peak at 
12007A. This would correspond to a rest-frame wavelength of 
1240A between the peaks of the N v doublet. The velocity off- 
set of the N v line from the systemic redshift Z[O mi] = 8.679, 


if indeed the 41238 line, is ~ +442 + 85 km s^!. This is 


consistent (within the error) to the offset we measure for the 
C ur] of ~ +310 + 150 km sl. If it is the 41243 line, the 
velocity offset is  —521 + 85 km s^!. 


4.4. Other Potential Broad Emission Lines 


As noted above, broad-line (Type 1) AGN typically ex- 
hibit broad permitted lines and narrow forbidden lines, in- 
terpreted as high-velocity dense gas in a “broad-line region" 
near the black hole (and dominated by its gravitational influ- 
ence) and lower-velocity gas in a “narrow-line region" that is 
more distant from the black hole and dominated by the galaxy 
kinematics (e.g. Netzer 2015). Unfortunately, the Ho line is 
redshifted beyond the wavelength range of NIRSpec for our 
target (though it is observable by MIRI spectroscopy). The 
other lines in the Balmer series (Hy, Hô, etc) are too weak for 
strong constraints on the presence of broad components (see 
Figure 6). In AGN with a broad Hf line, broad features are 
also typically observed in the permitted UV lines Lya, C tiv, 
C ur], Mgu with similar or stronger fluxes to the broad Hf 
line (e.g., Vanden Berk et al. 2001). 

We detect a significant broad width in the observed Lya 
line of our target. The blue-side absorption of the Lya forest 
makes it difficult to measure the line width precisely, but we 
note that our asymmetric Gaussian fit (described in $ 3.3) 
finds a red-side width of FWHM=1360.6 + 479.9 km s7!: 
comparable to the best-fit width of the broad Hf component. 

The best-fit Gaussian to the C 1] line also indicates a broad 
line width of FWHM=643+157 km s^!. The observed C nı] 
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feature is not well-fit by two narrow Gaussians, with widths 
constrained to be the same as the [O 11]A5008 line for each 
of the doublet lines. The line width of the best-fit broad C ri] 
profile is somewhat narrower than would be expected for a 
stratified broad-line region, in which the higher-ionization 
C mı] gas orbits closer to the black hole and consequently has 
a broader velocity width than H. The line center of the 
best-fit broad C m] profile is redshifted by +310 km s^! with 
respect to the systemic redshift (determined from the narrow 
[O11]45008 line as described in § 4.1). The shift in line 
center and narrower-than-expected width may indicate a blue- 
shifted absorption component in the C nı] line, as commonly 
observed in broad absorption-line quasars (e.g., Gibson et al. 
2009). 

We additionally see a significant broad component in the 
best-fit profile for the N1v]A1486 line. This feature is not 
typically observed in luminous quasars at z < 4, but is de- 
tected as a narrow feature in the NIRSpec prism spectrum 
of the z = 10.6 galaxy GN-z11 (Bunker et al. 2023). Its 
presence in GN-z11 has been interpreted as evidence for 
high-density and nitrogen-enhanced gas (Cameron et al. 2023; 
Senchyna et al. 2023), both properties consistent with expec- 
tations for broad-line regions around rapidly-accreting and 
low-mass AGN (Hamann & Ferland 1999; Matsuoka et al. 
2017). 

On the other hand, we do not observe broad components 
in the Civ and Mgu emission lines. The absence of Civ 
can be explained by the low signal-to-noise in its region of 
the NIRSpec spectrum, and the large flux uncertainty at that 
wavelength range allows for an undetected broad C rv line that 
is of comparable brightness to the broad H£ feature. The Mg 11 
feature is more puzzling: the observed Mg u profile is best-fit 
by narrow Gaussians for each of the doublet lines, consistent 
with galaxy (or narrow-line region) emission rather than the 
expected broad feature of an AGN. It is possible that this line is 
intrinsically broad but is affected by intervening absorption by 
the Mg u doublet: this requires a particular blueshift velocity 
of the absorption component but is not uncommon in quasars. 
Alternatively, the AGN may be attenuated by dust such that 
its UV emission lines are weaker than the Balmer lines, as 
implied by the best-fit SED model, although this conflicts with 
the detection of broad Lya, N 1v], and C ni]. Deeper spectra of 
this source would more effectively test for the presence of the 
broad UV lines expected for AGN, especially in investigating 
the hypothesis for blueshifted absorption affecting the broad 
C ni] line and the Mg i1 line. 


4.5. Emission-Line Flux Ratio Measurements 


Despite potential lingering issues in the absolute flux cal- 
ibration of the NIRSpec instrument, we find that the rela- 
tive flux calibration is consistent at least for pairs of emis- 
sion lines near one another in wavelength. We measure 


Table 4. Emission Line Ratios for CEERS. 1019 


Emission Line Ratio Measured Value Convention 


Oxygen 
[O 11]500$ / [O 11]4960 3.03 + 0.13 
[O 11]4363 / [O 11]500844960 0.02 + 0.01 RO3 
[O III]5008 / [O 11]372743729 13.46 + 1.18 032 
O/H. 
[O 11]5008+4960+10 111372743729 / H6 14.67 + 2.68 R23 
[O 11]372743729 / Hg 0.78 + 0.16 
[O 11]372743729 / Hy 1.13 + 1.14 
[O 11]4363 / Hy 0.45 + 0.16 
[O u1]500g / HB 10.44 + 1.91 R3 
Ne 
[Ne u1]3869 / [O 111372743729 1.18 + 0.14 
[Ne u1]3869 / Hy 1.34 + 0.22 
Balmer Series 
Hy / HE 0.69 + 0.16 

0.36 + 0.06 Total HB 
Hf /Hy 1.46 + 0.34 

2.82 + 0.49 Total HB 
Hf / Hô 3.33 + 1.24 

6.42 + 2.18 Total HB 
H8/He 2.43 € 0.76 

4.69 + 1.29 Total HB 


Nore—Common emission line flux ratios as measured from 
JWST/NIRSpec for CEERS_1019. All of the ratios that include 
H£ use the narrow component flux from the Dual Component Gaus- 
sian fit (§3.4.3) as reported in column 6 of Table 4 unless otherwise 
noted as using the Total H£ which is the value in column 4. 


[O 111]5008/[O i11]4960. = 3.03 + 0.13 which is consistent with 
the atomic physics calculation within errors (Storey & Zeip- 
pen 2000). We are confident in using the ratios of emission 
lines close in wavelength. Still, we acknowledge that ad- 
ditional flux calibration and improvements of the NIRSpec 
instrument are needed to trust widely separated line ratios 
fully. We report them in this paper given the above caveat. 

We list the measured line ratios for this source in Table 4. 
All reported line ratios are measured without using the Broad- 
line fit to the Hf line, and only the narrow component from 
the Dual component fit (See Figure 6). For any of the line 
ratios reported where one line is not detected at a significant 
level, we report the 1c upper limit for the line and mark it with 
an asterisk in the tables of Line Fluxes and Ratios (Tables 3) 
and 4). 
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Figure 7. The images show 5" x 5" cutouts around this source in the CANDELS HST/ACS optical bands, CEERS NIRCam near-infrared 
bands, and CEERS MIRI bands. The source is a clear dropout in the optical and is well detected at 1-8um. The large inset plot shows the 
photometric SED from HST and JWST. We also show Spitzer/IRAC measurements from Finkelstein et al. (2022a) with both TPHOT (light 
gray) and Galfit (dark gray); this source was highly blended in IRAC, and the Galfit measurements appear closer to the NIRCam measurements. 
The lines are best-fit models from Prospector (gray) and BAGPIPES (purple). The SED is dominated by stellar emission ($5), and these stellar 
models infer a massive (log M/Mo ~ 9.5) and highly star-forming (log sSFR ~ —7.9 Gyr-!) stellar population. The small inset panel shows 
constraints on the photometric redshift before (green) and after (red and yellow) the inclusion of JWST data. All three are consistent with the 
spectroscopic redshift, with those including JWST data placing much tighter constraints. 


5. CONSTRAINTS FROM CONTINUUM EMISSION 


We explore fitting the photometric SED of this object with 
stellar and AGN emissions simultaneously to explore which 
dominates the observed continuum emission. We first use 
the FAST v1.1 (Kriek et al. 2009; Aird et al. 2018) SED 
fiting code, including a star-forming galaxy and an AGN 
component as described in Kocevski et al. (2023). With this 
code, we try fits using two sets of AGN templates. The first 
uses eight empirically determined AGN templates. These in- 
clude five AGN-dominated templates from the Polletta et al. 
(2007) SWIRE template library (namely, the Torus, TQSO1, 
BQSOI, QSOI and QSO2 templates) and three composite 
SEDs of X-ray-selected AGN with absorption column densi- 
ties of Ny = 1022-23, 1023-24 and 1024-25 cm~? from Silva 
et al. (2004). See Appendix A of Aird et al. (2018) for addi- 
tional details. We also try a second fit with a low-metallicity 


(Z € 0.4 Zo) AGN component from CLoupy modeling, cho- 
sen to mimic the SED of a radio-quiet AGN (see $6.5 for 
more details). We measure the ratio of light from the best- 
fitting AGN model to the total stellar+AGN model in two 
wavelength windows: a rest-UV window at 0.15-0.25 um, 
and a rest-optical window at 0.51-0.6um (both designed to 
avoid strong emission lines). In both of these fits, the pho- 
tometry constrains the AGN to be subdominant to the stellar 
emission. With the fiducial AGN template, the stellar emis- 
sion comprised 99.5% of the UV flux and 85% of the optical 
flux. With the low-metallicity AGN model, stellar emission 
comprises 82-83% of both the UV and optical flux. 

We try an alternative fit using the Cigale code (Boquien 
et al. 2019; Yang et al. 2020, 2022). In the fit, we adopt a 
modified Schartmann et al. (2005) AGN accretion disk model. 
We allow the deviation from the default optical spectral slope 
(OaGn parameter) varying from —1 to 1 and polar-dust extinc- 
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Table 5. Physical Properties of CEERS_1019 


SFH Prior Muy,obs UV Slope 6 log M,  logSFRso log sSFRs5q Ay log Z 
[mag] [Meyr!] [ye] [mag] ^ [Ze] 
(1) (2) (3) (5) (6) (7) (8) 
Continuous -22.07:00 1.774045  94*05 — 1.5403 -719:02 0402  -0.3:05 
Bursty — -22.00'00 — —1,80'009 9.5302 — 1,5304  — —79*02 04705  -0.6396 


NorE— Physical properties of CEERS 1019 as derived through SED fitting to the HST and JWST photometry at the spectroscopic redshift. 
SED fitting with Prospector (Johnson et al. 2021) using both a continuity and bursty SFH prior as described in Tacchella et al. (2022). (1) 
Star-formation history prior used in fit. (2) Absolute UV magnitude at rest-frame 1500À. (3) UV Spectral slope as measured from the model 
spectra. (4) Stellar mass. (5) SFR averaged over the past 50 Myr (6) Specific SFR averaged over past 50 Myr (7) Dust attenuation at 5500, from 
SED models. (8) Stellar metallicity from SED models to the HST photometry. (4—8) values are calculated from SED fitting as in Tacchella et al. 


(2022). 


tion [E(B — V)pp parameter] varying from 0 to 1 (see 84 of 
Yang et al. 2022) for details. Measuring the results in the 
same two wavelength windows as above shows that the stellar 
emission comprises 70% and 59% of the continuum emission 
in the rest-UV and rest-optical, respectively. 

These results imply that the observed continuum emission 
is dominated by stellar light. We thus proceed to perform 
Bayesian SED modeling with stellar-only models to explore 
what constraints can be placed on the stellar population of the 
host galaxy. Our fiducial results come from the Prospector 
SED fitting code (Johnson et al. 2021), following Tacchella 
et al. (2022) using both a continuity and bursty prior on the 
star-formation history (SFH). Both SFH priors lead to similar 
posterior distributions, implying that the data is informative 
and the impact of the prior is minimal. We find that this 
galaxy has a stellar mass of M, ~ 10°>*°3 Mo and is actively 
forming stars (SSFR ~ 10772493 yr-!), doubling its mass 
every ~ 100 Myr. We find a mass-weighted age of tso = 
34*)P Myr. We list the physical properties of this galaxy 
from our SED Fitting in Table 5 and show the model fits in 
Figure 7. 

We also fit these data with BAGPIPES (Carnall et al. 
2018). We generally followed the procedures in Papovich 
et al. (2022), but we have now used BPASS v2.2.1 stellar pop- 
ulation models (Eldridge et al. 2017), a star-formation history 
represented as a Gaussian Mixture Model (Iyer et al. 2019), 
and fit over a range of ionization parameters log U = [—4, —1] 
to model the nebular emission. From these fits, the galaxy has 
a stellar mass of log M./ Mo = 9.3 + 0.1, consistent with the 
fits from Prospector, though a slightly higher log(sSFR) of 
-7.5 x 0.2 yr^!. The systematic uncertainties on stellar mass 
and SFR here are «0.2 dex. Inferred properties, such as stel- 
lar mass, SFR and dust attenuation, are generally consistent 
and within the systematic uncertainties of the Prospector fits, 
but require a slightly higher SFR to account for differences in 
modeling assumptions. We show the best-fitting Prospector 
and BAGPIPES models in Figure 7, alongside images of this 
galaxy in all observed filters and photometric redshift results 
from EAZY. 


5.1. UV Slope (B) from Photometry 


To measure the UV spectral slope, 5, for this source, we 
fit a power-law function (fy œ AP; Calzetti et al. 1994) to 
the observed photometry. We only include filters within the 
rest-frame 1500 — 3000 À range to avoid contamination from 
the Lya@ break and strong emission lines. For this source, 
these filters are F160W, F150W, F200W, and F277W. Us- 
ing the EMCEE software (Foreman-Mackey et al. 2013), we 
measure the posterior distribution on f. CEERS 1019 has 
B= —1.76+9:12, where the uncertainty is taken as the 68% 
central width from the posterior. This is consistent at the 
1c level with the value inferred from HST photometry from 
Tacchella et al. (2022) of B = -1.6170.18. It is also consis- 
tent with the value of 8 measured from the Prospector model 
spectra, listed in Table 5. This value of 8 is comparable to 
the UV colors from other similarly bright sources (Tacchella 
et al. 2022), consistent with the low but non-negligible dust 
attenuation we find from our SED fit (Ay = 0.4 + 0.2). 


5.2. Morphology 


We investigate the morphology of this source using the fit- 
ting codes Galfit (Peng et al. 2002, 2010) and statmorph 
(Rodriguez-Gomez et al. 2019). We use the Galfit least- 
squares fitting algorithm to fit the galaxy's surface brightness 
profile. We use a 100x100 pixel cutout of the F200W science 
image as input, the corresponding error array (the ‘ERR’ ex- 
tension) as the input sigma image, and the empirically derived 
PSFs. As an initial guess, we use the source location, magni- 
tude, size, position angle, and axis ratios from the SE catalog. 
The Sérsic index is allowed to vary between 0.01 and 8, the 
magnitude of the galaxy between 0 and 45, the effective radius 
(Re) between 0.3 and 100 pixels, and the axis ratio between 
0.0001 and 1. We also allow Galfit to oversample the PSF 
by a factor 9. We then visually inspected the best-fit model 
and image residual for each source to ensure that the fits were 
reasonable and that minimal flux remained in the residual. 

We find that the source is poorly fit by a single Sérsic 
component. Visually, the source is extended asymmetrically, 
with three distinct components (see Fig. 8). The optimum fit 


16 LARSON ET AL. 


F356W 


Figure 8. JWST NIRCam 1.5" x 1.5" cutouts of CEERS 1019 in four different filters (F200W, F277W, F356W, and F444W) and an RGB 
color composite (with blue=F115W+F150W+F200W, green=F277W+F356W+F410M, and red=F444W) with each filter at its native resolution, 
highlighting the substructure visible at shorter wavelengths. The positions of the NIRSpec MOS shutters are overlaid. This source has a bright 
central component centered in the shutter and two extended components as discussed in $5.2. 


is obtained when the central region is fit with both a point 
source and a Sérsic component, and the two other regions to 
the west and the northeast are each fit with their own Sérsic 
component. The requirement of a point source component 
for the fit is consistent with the source having an AGN. We 
similarly performed fits on the F356W and F444W images and 
again find that single Sérsic fits do not work because of the 
multiple components, even with the lower resolution at these 
longer wavelengths. The asymmetric nature of the source and 
the presence of the three separate components are consistent 
with the galaxy being involved in a major merger. We use 
statmorph to measure the concentration parameter and size 
as a function of wavelength on PSF-matched cutouts and find 
that C=2.82 and R,26.49 pixels (0.91 kpc) for F356W and 
C=2.93 and Re=5.64 pixels (0.79 kpc) for F444W, suggesting 
that the F444W light is more concentrated than the F356W or 
F200W emission, supporting the visual impression in Figure 
2. At the redshift of this galaxy, the F444W emission is 
dominated by the [O mm] emission line (see Fig.7) 


6. CONSTRAINTS FROM NEBULAR LINE EMISSION 
6.1. Black Hole Mass 


In this section, we estimate the virial mass of the central 
Black Hole (BH), assuming that the observed broad Hf emis- 
sion traces the kinematics of gas in the broad-line region. The 
measurement of BH masses with single-epoch spectra can be 
made using the width of the broad Hf emission line and the 
rest-frame 5100 À continuum luminosity, L51900, which has 
been shown to correlate with the distance to the broad-line 
region using reverberation mapping (e.g., Kaspi et al. 2000; 
Cackett et al. 2021). However, this assumes the rest-frame 
5100 À continuum luminosity is dominated by light from the 
AGN, which is likely not the case here (see $5.1). As a result, 
we instead use Equation 10 from Greene & Ho (2005), which 
employs only the H£ line luminosity, Lyg, and width: 


Lug | FWHMmgg ) M 
ion 


Mpu = 2.4 x 10° 
= | 10? erg s7! 10? km s7! 


(5) 


Table 6. CEERS 1019 Black Hole Properties 


log Mgy 6.95 + 0.37 [Mo] (1) 
FWHMyg, Broad — 1196.26 + 349.07 [kms7!] (2) 
log Lug, Broad 42.5 € 02 [ergs-!] (3) 
log Lgol 45.1 + 0.2 [ergs] (4) 
AEdd 1.3 € 0.5 (5) 


NorE— (1) Black Hole Mass (86.1). (2) FWHM of the broad Hf 
line ($4.2). (3) Luminosity of the broad Hf line. (4) Bolometric 
luminosity. (5) Eddington ratio (Lpoi/ Leda). 


This equation is based on the formula of Kaspi et al. (2000), 
but with Lyg substituted for 75,09 using the empirical corre- 
lation between Balmer emission-line luminosities and 5100 
reported in Greene & Ho (2005). 

Using the luminosity and width of the broad Hf compo- 
nent, we derive a BH mass of log(Mgy/Mo) = 6.95 + 0.37. 
This mass is 1-2 dex lower than existing samples of luminous 
quasars with BH mass estimates at z > 5 and more compa- 
rable to the low-luminosity AGN reported in Kocevski et al. 
(2023), found in these same CEERS NIRSpec data. Our mea- 
sured mass implies that the AGN in CEERS 1019 is powered 
by the least massive BH known in the Universe during the 
epoch of reionization. 

To determine the accretion rate onto the BH relative to the 
Eddington limit, we derive a bolometric luminosity, Lgo from 
the broad H line luminosity and compare it to the Eddington 
luminosity, Lgaa, for our measured BH mass. Assuming an in- 
strinsic broad-line Ha /H£ ratio of 3.06 (Dong et al. 2008) and 
a bolometric correction of Lg) = 130 x Lye (Richards et al. 
2006; Stern & Laor 2012), we obtain Lgo = 1.4 + 0.6 x 10° 
erg sl. This results in an Eddington ratio, Agga = Lpoi/ Leda, 
of 1.3 + 0.5. This suggests that the BH is undergoing a rapid 
growth phase and is accreting at approximately its Eddington 
limit. 
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6.2. Electron Temperature and T,-Based Metallicity 


The [O m1]4364/[O u1]500844960 ratio can be used to mea- 
sure the electron temperature (Te) of the galaxy's interstellar 
medium (ISM). As these lines are all collisionally excited, 
but since the [O rir]43g4 line de-excites at a higher level, a 
higher [O 111]4364 emission relative to [O 111]5008+4960 indicates 
higher-energy electrons are responsible for the excitation. 
This Te has been found to correlate with the ISM metallicity, 
providing a way to 'directly* measure the metallicity of the 
source with these lines (i.e., Kewley et al. 2019b). 

Our source has [O m1]4364/[O 1n1]4960+5008 = 0.022 + 0.007. 
We use the measured [O mi]4364/[O 111] 496045008 ratio to esti- 
mate 7, using Equation 4 of Nicholls et al. (2020): 


3.5363 + 7.2939x 


logs 6 
o810(Te) = T 5000 1.6298x 0.12213 — 000743. © 
where 
[O 111] 4364 | 
x =log,, | 4 — | = -1.657 + 0.003 
Fin (a 111]4960+5008 


This yields logįoọ (Te) = 4.270 + 0.566, and T, = 18630.755 + 
3.682K for CEERS_1019. 

By measuring 7, in this way, we are only sensitive to the 
portion of the ISM emitting the [O m] lines, which may not 
represent the entire galaxy. Any significant gradients in den- 
sity and/or ionization may lead to a mixture of different re- 
gions in the galaxy, and the T, from [O m1] is only sensitive to 
the high-ionization regions. 

We then use Equation 1 from Pérez-Montero et al. (2021) 
to estimate the metallicity from T, as shown below: 


12 + log(O/H) = 9.72 — 1.70te + 0.3222 (7) 


where t, = T, in units of 104K and thus our ‘direct’ metallicity 
measurement for this source is: 12+log(O/H) = 7.664 + 0.508, 
or 0.095*0205 Zo (Zo=8.69, Asplund et al. 2021). 

This method is used for extreme emission line galaxies 
(EELGs), galaxies with extreme emission lines dominated by 
stellar photoionization, so it may not be strictly accurate for 
objects whose photoionization is dominated by an AGN. The 
above equation from Pérez-Montero et al. (2021) is adopted 
from Amorin et al. (2015), which derived this relation based 
upon calibrations using local EELGs, after discarding AGN 
candidates based upon BPT (Baldwin et al. 1981) diagram 
measurements. 


6.3. Electron Density from [O 11] Doublet 


The ratio between the two lines in the [O 11]372743729 doublet 
is often used to determine the electron density (ne) in star- 
forming regions at temperatures T~ 10,000 = 20, 000K. This 
is because the excitation energy between the lines is on the or- 
der ofthe thermal electron energy. Thus the relative excitation 


rates depend only upon collision strengths (Osterbrock 1989). 
While the [O 11]372743729 doublet is blended at the resolution 
of the NIRSpec data (R~1000 in the M gratings), the separa- 
tion between the two line peaks at this redshift is large enough 
to distinguish. The Doublet fit to this line ($3.4.2) gives flux 
values for both emission lines such that we can estimate a 
ratio between the two of [O 11]3729 / [O 11]3727 = 0.639 € 0.255 
with which to infer an electron density, n, = 1.9 € 0.2 x 10? 
cm^? (Osterbrock 1989). Alternatively, measurements using 
the Cr] doublet are sensitive to larger ne. We only de- 
tect the C imi]19os line from this doublet which could imply 
that CEERS. 1019 has very high density, ~ ne > 10* cm? 
(Keenan et al. 1992). 


6.4. A, from Balmer Decrement 


The ratios of Balmer lines can be used to measure nebular 
attenuation. In particular, gas with Case B recombination and 
a temperature of T, = 20, 000K (as calculated for our object in 
§ 6.2) will have intrinsic Balmer ratios of H/Hy = 2.110 and 
Hjf/Hó = 3.817 (Osterbrock 1989). We can estimate the dust 
attenuation E(B-V) with an assumed attenuation curve (here, 
we use Calzetti 2001) by comparing our observed line ratios 
to these intrinsic values. For the observed line ratios we use 
the total line flux, rather than the narrow-only flux, because 
the Hy and Ho lines do not have sufficient signal-to-noise to 
separate broad and narrow components. This means that our 
estimated nebular attenuation may be a mix of attenuation 
affecting the AGN and the narrow lines. 

From the observed line ratios of HB/Hy = 2.82 + 0.49 and 
Hf/Hoó = 6.42 + 2.18, we find nebular attenuation of: 


* E(B-V)yp/Hy = 0.64 + 0.41 
* E(B - V)gg;gs = 0.83 + 0.54 


The observed line ratios are statistically consistent with the 
intrinsic values such that the implied E(B — V) values are 
similarly consistent with little to no dust. But the large error 
bars also allow for a broad range of nebular attenuation. In 
§ 5.1, the Prospector SED fit implied Ay =0.4 + 0.2, which 
equates to E(B-V) = 0.1 + 0.05. This is consistent within 
the large error bars of the Balmer-decrement estimate for dust 
attenuation. 


6.5. Ionization Parameter? 


To place this source in the context of other star-forming 
galaxies and AGN identified in this early epoch, we investigate 
the “OHNO” line ratio diagnostic diagram. This diagram, 
comparing [O m]s008/H versus [Ne 111]3g70/[O 11]3727+3729, 


^ Defined as the ratio of the number density of ionizing photons to 
the gas density, U = ny /ng. 
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Figure 9. The “OHNO” diagram, the line ratio diagnostic using [O m]5008/H£ versus [Ne r1]3870/[O 11]37274.3729. with the black curve showing 
the boundary between star-forming and AGN (as defined in Backhaus et al. 2022, with the shaded region showing the star-forming area). The 
orange stars are from this work, representing the narrow component fit (big star) and single Gaussian fit (to the full profile; small star) of the HB 
emission line in our galaxy. The black polygons represent other high-redshift (z > 5) galaxies with JWST/NIRSpec observations (Trump et al. 
2022; Tang et al. 2023; Kocevski et al. 2023). Underplotted on both panels are photoionization models, with the left panel showing star-forming 
(blue) and AGN models (red) using the CLoupv v17 code (Ferland et al. 2017) and the right panel showing star-forming models (green) from 
Kewley et al. (20192) using the MAPPINGS V code (Sutherland et al. 2018). Each line of points represents a specific set of ionization parameters 
of logioU = [-2.5,-2.1,-1.5] for the CLoupv models and log;9U = [—2.5,—1.5] for the MAPPINGS models, with increasing ionization towards 
darker colors (or, towards the top right of each panel). For the MAPPINGS models, the lines are further separated by isobaric pressure (denoted 
by line type). The metallicities of the models are represented by the size of the points. 


has been used at lower redshifts to identify ionization pow- 
ered purely by star formation from that of an AGN (e.g., 
Backhaus et al. 2022; Cleri et al. 2022, 2023). Recent studies 
of high-redshift galaxies (z » 5) using JWST spectroscopy 
have found most if not all of their sources occupying AGN 
region of the diagram (i.e., above the theoretical line dividing 
star-formation from AGN, calibrated via low-redshift obser- 
vations) — results which seem to suggest strong ionization 
from AGN and/or very metal-poor star-forming H II regions 
at these high redshifts (e.g., Trump et al. 2022; Tang et al. 
2023; Übler et al. 2023; Kocevski et al. 2023). 

Such is the case for CEERS_ 1019, with high measurements 
of the ionization indicators [O un1]5008/H£ = R3 = 10.44 + 1.91 
and [O 111]4960+5008/[O 11]372743729 = O32 = 13.46 + 1.18 (see 
Table 4). Figure 9 shows our galaxy on the “OHNO” dia- 
gram (orange stars), using the narrow component of the H8 
flux (bigger star; also the light green line in Figure 6; the 
smaller star shows the ratio with the full Hf flux). Additional 
sources from recent JWST studies focused on the reioniza- 
tion era (Trump et al. 2022; Tang et al. 2023; Kocevski et al. 
2023) are included as small gray polygons. The black curve 
in both panels denotes the boundary between star-forming 
and AGN as defined by Backhaus et al. 2022, with the associ- 
ated grey-shaded region showing the z ~ 1 “star forming” re- 
gion. Under-plotted are photoionization models from both the 


CLoury v17 (Ferland et al. 2017, left panel) and MAPPINGS 
V (Sutherland et al. 2018, right panel) codes, with parameters 
chosen to showcase a range of possible galaxy properties at 
this redshift, including elemental abundances, stellar popula- 
tions, and stellar and gas-phase metallicities. In this figure, 
each line in the left panel represents CLoupy models with a 
specific ionization parameter (increasing logio U shown by 
the darker colors), where the size of points in each line repre- 
sents the specific gas-phase metallicity of said model. Simi- 
larly, each line in the right panel represents the MAPPINGS 
models for a specific ionization parameter (similar to the left 
panel) and isobaric pressure (denoted by the line types). 

The CLoupy models shown in Figure 9 assume a plane- 
parallel geometry with a nebular electron density of 10? cm? 
(inferred from the flux ratio of the [O 11]372743729 doublet, see 
$6.3), scaled Solar elemental abundances, and cover ioniza- 
tion parameters from logio U = [-2.5, -2.1, -1.5]. The 
star-forming models use the Binary Population and Spectral 
Synthesis v2.0 (BPASS, Eldridge et al. 2017) fiducial binary 
stellar population models with an IMF that extends to 300 
Mo (described by a slope of a = —1.30 from 0.1-0.5 Mo 
and œ = —2.35 from 0.5—300 Mg) and continuous star for- 
mation of 1 Mọ yr^!. For these models, we fix the stellar 
and gas-phase metallicities, covering Z = 0.1-0.5 Zo. The 
AGN models use the table agn model in CLoupy, which 


z=8.679 AGN witH CEERS 19 


approximates a "typical" radio-quiet AGN, covering the same 
range of ionization parameters as the star-forming models and 
spanning nebular metallicities of 0.05-0.5 Zo. 

The MAPPINGS models shown in the same figure are the 
“Pressure Models” from Kewley et al. (2019a), assuming a 
plane parallel geometry with a range of pressure (logig(P/k) 
= 7-9 cm^?), ionization parameter (logjgU = [-2.5,-1.5]; 
also referred to as log19Q = 8-9), and metallicity (Z = 0.05—1 
Zo). These star-forming models use the Starburst99 (Lei- 
therer et al. 2014) stellar population models with a Salpeter 
IMF (Salpeter 1955) extending to 100 Mo. The stellar and 
gas-phase metallicities follow a prescription such that at lower 
metallicities, the models have increasingly enhanced «o abun- 
dances (see Nicholls et al. 2017 for more details). We have 
included these models in addition to the CLoupy models 
(which have Solar abundances scaled to chosen metallicity) 
to showcase the range of potential galaxy properties for our 
high-redshift galaxy. 

Similar to that found in other spectroscopic studies in this 
epoch (e.g., Trump etal. 2022; Tang et al. 2023; Kocevski et al. 
2023), the line ratio when using the single component fit to the 
Hf line for CEERS 1019 (the smaller star) lies in a region on 
this diagnostic which is difficult to differentiate between AGN 
and metal-poor high ionization star-forming H II regions. 
This is not unexpected, as galaxies at higher redshifts have 
been shown to generally have higher ionization and lower 
metallicities in comparison to those at lower redshifts (e.g., 
Shapley et al. 2003; Erb et al. 2010; Kewley et al. 2019a; 
Backhaus et al. 2022; Papovich et al. 2022; Sanders et al. 
2023). This suggests that at these high redshifts, successfully 
separating star-forming and AGN sources using these strong 
line diagnostics can be challenging — evidenced by the lower 
metallicity AGN models overlapping with the star-forming 
models in the left panel of Figure 9. Further discussion about 
the general utility of such emission line diagnostics at high- 
redshift is discussed in $7.3. 

However, when focusing on the fit to only the narrow com- 
ponent of the Hf line for CEERS. 1019 (the bigger star), there 
is a clear distinction from the rest of the high-redshift sources 
shown in this figure. Of the sources shown, our galaxy is 
the only one with enough S/N in Hf to a) see a broadening 
of the line profile and b) clearly measure both a broad and 
narrow component for the line.? Follow-up spectroscopy of 
these and other sources that fall in this strong line regime at 
high-redshift may shed light on more "hidden" AGN in the 
early Universe (Kocevski et al. 2023). 

From Figure 9, CEERS 1019 covers a similar location to 
the CLoupy star-forming models with logi9U ~ —1.9 to —1.5 


5 In Kocevski et al. (2023), while they could resolve both components for the 
Ha emission in the two AGN in their sample, the HG emission for both 
sources had S/N too low to enable this identification. 


(for the single component fit to Hf, smaller star) and the 
AGN models with log;9U ~ —2.1 (for both the narrow and 
single component fits to H£; both stars). Similarly, the MAP- 
PINGS models suggest agreement with logigU ~ —2.5 to 
—1.5; however, it is also dependent upon the pressure as- 
sumed. These results add to the expectation that this galaxy 
has high ionization powering these strong line ratios. Indeed, 
this agrees well with the strong O32 line ratio measured for 
this galaxy which is often indicative of highly ionized, metal- 
poor gas (e.g., Schaerer et al. 2022; Williams et al. 2022) and 
a regime that could suggest Lyman continuum leakage (e.g., 
Izotov et al. 2018). As a comparison to the ionization pa- 
rameters gleaned from Figure 9, we estimate this value using 
two relations from the literature. First, using the theoretical 
relation between O32, 12+log(O/H), and ionization param- 
eter described in Kobulnicky & Kewley (2004), we derive 
logio ^ —2.3. Finally, using the measured relation between 
O32 and the ionization parameter explored in Papovich et al. 
(2022), we derive logio ~ —2.0. These results are all rela- 
tively consistent with one another and further highlight the 
strong ionizing nature of this source. 


7. DISCUSSION 


The "chicken or egg" origin of BH seeds in the first galax- 
ies remains unsolved. Theoretical predictions suggest a mix 
of "light" (^10? Mo) Population III stellar remnants and/or 
"heavy" (^10? Mo) seeds formed via the direct collapse of 
primordial gas, as summarized in recent reviews by Smith & 
Bromm (2019), Inayoshi et al. (2020) and Fan et al. (2022), 
but the relative mix of each seed type remains unconstrained 
by observations. Regardless of the seeding mechanisms, our 
Universe is capable of forming extremely massive SMBHs 
very early, with the highest-redshift massive quasar at the 
time of this writing being z — 7.642 (Wang et al. 2021), exist- 
ing $700 Myr after the Big Bang. Observations with JWST, 
especially MIRI observations of obscured growth (e.g., G. 
Yang et al. in prep) and NIRSpec observations of broad (z 
1000 km s^!) and/or high-ionization (e.g., N V, He II, [Ne V], 
C IV) emission lines can now allow the first real census of 
AGN in low-mass (M, « 10!? Mo) hosts at high redshift, 
providing strong constraints on the initial black hole seed 
distribution. Direct observations of early SMBHs will not 
only constrain mechanisms for early black hole growth (e.g. 
Ricarte & Natarajan 2018a,b), but can also provide further 
insight into the role ionizing photons from AGN played in the 
reionization of the IGM (e.g., Finkelstein et al. 2019; Gial- 
longo et al. 2019; Dayal et al. 2020; Grazian et al. 2020; Yung 
et al. 2021; Grazian et al. 2022). 


7.1. Constraints on the Formation of this z=8.68 AGN 


Discovering an actively accreting SMBH at such a high 
redshift further constrains the formation time for such ob- 
jects. Here we consider plausible formation mechanisms 
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Figure 10. Mass of black holes versus redshift. The gold star represents the result of this work. The green points denote published z = 6—7.64 
quasars, taken from the compilation of Inayoshi et al. (2020) augmented by more recent z > 7 quasars from Fan et al. (2022). The different 
lines show potential tracks of black hole growth. The purple lines show a 100 Mo Population III stellar remnant seed forming at z ~ 30, with 
growth beginning after a 100 Myr delay due to progenitor gas heating. The dashed purple line shows constant Eddington growth, which is both 
likely un-physical and cannot reproduce our observed object mass. The solid line shows one example of how periods of super-Eddington growth 
(10x Eddington for 10 Myr) separated by periods of sub-Eddington growth (0.1x Eddington for 50 Myr) could plausibly create the observed 
source. The red lines show potential formation mechanisms from a small (3x 104 Mo) DCBH forming at z = 15 (dotted), and a large (10° Mo) 
DCBH forming at z = 10.5 (dashed), each growing at the Eddington rate. Either of the mechanisms that could explain our observed black hole 
source, DCBH+Eddington or stellar+Super-Eddington, and are somewhat exotic scenarios, pushing standard assumptions. Probing SMBHs out 
to higher redshifts and lower masses will clarify the formation mechanisms of these objects. 


for this source (e.g., small versus large seeds), as well as 
discuss what its plausible descendants are. While the in- 
ferred mass of this BH is not much larger than the proposed 
range for DCBH seeds, it is unlikely that we are witnessing 
a DCBH soon after formation. The constraints we place on 
both the stellar mass and gas-phase metallicity of the host 
galaxy are significantly higher than expected (DCBH forma- 
tion requires near-primordial conditions), implying that this 
object is a "standard" (albeit very distant) AGN accreting at 
~the Eddington rate, with a continuum SED dominated by 
stellar emission (in line with the stellar emission dominating 
the continuum SED, §5), similar to the scenarios recently 
explored by Volonteri et al. (2023). 

It is therefore interesting to explore how this AGN was 
seeded: from a low-mass (~ 10-100 Mo) stellar seed or a high- 
mass (^ 1046 Mo) DCBH seed, summarized in Figure 10. 
The purple curves show idealized black hole mass tracks 
assuming a 100 Mo stellar seed. Such a seed could plausibly 
form at z ~ 30 but would be unable to accrete for ~100 Myr 
due to radiative heating of the gas from the stellar progenitor 
(e.g., Johnson & Bromm 2007; Jeon et al. 2014). Assuming 


such a 100 Myr delay, the onset of black hole growth would 
begin by z ~ 18.5. As shown by the dashed purple line, 
Eddington-limited accretion in this scenario would be unable 
to reach the inferred black hole mass by the redshift measured 
for this object. The enhanced stellar feedback environment 
the host galaxy must have experienced in the previous ~100 
Myr further complicates this. Building up to the observed 
Z z 0.01Ze requires a large number of supernova explosions 
in a short period of time. This would lead to violent, turbulent 
mixing and heating of the gas, making accretion likely very 
inefficient (significantly sub-Eddington). 

However, there may also be periods of super-Eddington 
(“catch-up”) accretion; these would occur in short episodes, 
with minimal sub-Eddington accretion between them; such 
short bursts of black hole growth could plausibly build up 
the «107 Mo SMBH we observe in an otherwise stellar- 
emission dominated galaxy (e.g., Volonteri & Rees 2005; 
Madau et al. 2014; Inayoshi et al. 2016). For example, Pez- 
zulli et al. (2016), when examining the evolution of merger 
tree simulations, find that super-Eddington accretion modes 
(e.g., Haiman 2004; Silk 2005; Polletta 2008) can be a signif- 
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icant component of SMBH growth in gas-rich environments, 
where up to 75 percent of the SMBH growth can be accounted 
for by periods of super-Eddington accretion, with intermittent 
phases of disruption resulting from rapid depletion / replen- 
ishment of the bulge gas reservoir out of which the BHs 
accrete. 

In Figure 10, we highlight one potential version of this 
scenario, with the purple solid line showing an object which 
starts with a seed mass of 100 Mo, which begins growing at 
z = 18.5 with episodic 10 Myr periods of super (10x) Ed- 
dington growth, followed by 50 Myr "breathing" periods with 
sub (0.1x) Eddington growth. Periodic super/sub-Eddington 
growth allows this hypothetical object to grow to ~ 107 Mo 
by z ~ 8.7 from a small seed, matching our observations. 
While plausible (in particular, the episodic growth with sup- 
pressed, sub-Eddington, accretion periods was predicted by 
the simulations of Jeon et al. 2012 and Massonneau et al. 
2023), this scenario is somewhat contrived, making a stellar 
seed somewhat unlikely. 

A DCBH origin may be more plausible. While the inferred 
metallicity implies the DCBH event occurred significantly 
prior to the observed epoch, this could provide the needed 
time for a massive seed ~10° Mo to grow by the needed fac- 
tor of ~10x without also having to invoke super-Eddington 
conditions. Given the likely strong stellar feedback environ- 
ment in the assembling galaxy, this DCBH scenario may be 
favored. We show two plausible DCBH tracks as the red lines 
in Figure 10. Both options, of a lower-mass 3 x 10* Mo 
DCBH forming at z = 15 or a higher mass 10° Mọ DCBH 
forming at z — 10.5, should they grow at the Eddington limit, 
could reproduce the observed mass of this object by z ~ 8.7. 

Lastly, extrapolating the inferred black hole growth tracks 
to lower redshifts provides some insight into the ultimate de- 
scendants of early black holes such as the one we observe 
here. As shown in Figure 10, this object is both too low-mass 
and observed too late to be the plausible progenitor of the 
z > 7 quasar population. However, should near-Eddington 
accretion be sustained for significant periods of time, objects 
like this one could plausibly evolve into the massive z ~ 6— 7 
SDSS quasar population. Further constraining both the pro- 
genitors and descendants of high-redshift SMBHs is possible 
as the population of z » 8 SMBHs is coming into view. Given 
the discovery of this source in a relatively small dataset, it is 
likely that further /WST spectroscopy, in particular, follow-up 
of bright sources discovered over wide areas with the upcom- 
ing Nancy Grace Roman Space Telescope, will yield such a 
population (e.g., Yung et al. 2023). 


7.2. Evolution in the Black Hole — Galaxy Mass 
Relationship 


Fig. 11 compares CEERS 1019 with other high-redshift 
(z ~ 6) AGNs (Izumi et al. 2021b; Kocevski et al. 2023) 


as well as the z = 0 Mpg — M. relationship. Because they 
are easier to detect with pre-JWST instruments, most high- 
redshift AGNs are high-luminosity. However, AGN have a 
well-known luminosity-dependent bias (i.e., the Lauer bias; 
Lauer et al. 2007), such that higher-luminosity AGNs tend to 
have over-massive black holes. 

JWST has improved our ability to select low-luminosity ac- 
creting black holes (e.g., Kocevski et al. 2023), which are 
expected to lie closer to the intrinsic (overall) Mgu — M. rela- 
tionship. Fig. 11 shows an estimate from the empirical TRIN- 
ıry model (Zhang et al. 2023) for how different bolometric 
luminosity thresholds bias the median Mpg — M. relationship. 
With an estimated bolometric luminosity of Lpo ~ 105 erg/s 
and a host stellar mass of M, ~ 10?2 Mo, TRINITY would 
suggest that the present source would be expected to have 
a AMpy ~ 0.7 dex offset from the median relation. The 
present source is located around the upper envelope of the 
Mggy — M, relation spanned by the AGN sample from Reines 
& Volonteri (2015) at z ~ 0. The present source also has a 
lower Mgy compared to the z = 5.55 AGN from Übler et al. 
(2023), which is ~ 10 times brighter. Qualitatively, this is 
consistent with the Lauer bias that brighter AGNs tend to be 
over-massive black holes compared to their host galaxies. 

If the present source lies on the median Mgy — M. for 
its luminosity, the Lauer bias estimate from TRINITY would 
suggest that the overall z = 8.7 median Mpy — M. would have 
Mpy = 1063 Mo at M. = 10?5 Mo, i.e., there would be 
no overall evolution from the z = 0 Mgu — M. relationship. 
If other sources are confirmed to have similarly consistent 
masses compared to the z = 0 Mpgy — M. relationship, it 
would cement a tight relationship between galaxies and black 
holes that extends into the reionization era, and would also 
provide a strong challenge to theoretical models that have 
predicted substantial evolution in the Mgy — M. relation for 
high-redshift and low-mass galaxies (see Habouzit et al. 2020, 
for a survey). 


7.3. The Efficacy of Line Ratio Diagrams to Identify AGN at 
Early Times 


Many of the canonical emission-line ratio diagnostics of 
AGN versus star formation using rest-frame optical emission 
lines are calibrated for low to moderate redshifts (z x 1, e.g., 
Baldwin et al. 1981; Veilleux & Osterbrock 1987; Trump 
et al. 2015; Backhaus et al. 2022; Juneau et al. 2011, 2014). 
As such, the divisions between the star-forming and AGN 
regions defined by these optical line ratio diagnostics are not 
necessarily valid at higher redshifts. 

Several works offer remedies to these low-redshift diagnos- 
tics by quantifying a redshift evolution to the division (e.g., 
Coil et al. 2015; Cleri et al. 2022), but these new divisions 
have not been extended to the EoR; in fact, these relations 
for these diagrams are shown to break down at high-z in sim- 
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Figure 11. The predicted z = 8.7 median M, — M. relation for quasars with different bolometric luminosity thresholds from the TRINITY model 
(solid curves; Zhang et al. 2023). The pink shaded region is the 1 — o spread around the median scaling relation for quasars (~ 0.55 dex), 
which includes the random scatter in observed M, when using virial estimates. This (log-)normal scatter is nearly luminosity-independent, so 
we only show it for the brightest quasars for clarity. The black solid line is the predicted Me — M. relation for all SMBHs at z = 8.7, and the 
black shaded region is the intrinsictobserved scatter around the intrinsic M, — M. relation. The green solid line is the M, — M, relation for 
AGNs brighter than log Lpoj [erg/s] 2 45.1, corresponding to the bolometric luminosity of the AGN in this work (gold star). For comparison, 
we also show the z = 0 relation from Greene et al. (2016) in the black dashed line, and the z = 0 AGN sample from Reines & Volonteri (2015) 
in black crosses. The following data are shown in symbols color-coded by bolometric AGN luminosities: 1) z ~ 6 quasars compiled by Izumi 
et al. (2021a) and Ubler et al. (2023) (filled circles); 2) The two z > 5 AGNs from Kocevski et al. (2023) (pentagons). 


ulations (Hirschmann et al. 2019; Hirschmann et al. 2022). 
The redshift-evolving diagnostics of these works employ a 
simple shift from the low-redshift divisions; however, there is 
insufficient knowledge of high-redshift AGN to validate this 
practice for galaxies in the EoR. 

For a different avenue to remediate these issues, other works 
have offered completely new emission-line ratio diagnostics 
with higher ionization emission-line ratios, e.g., in the optical 
with He II/H and [Ne V]/[Ne ui] (Katz et al. 2022; Cleri et al. 
2023), and in the UV with, e.g., C mi]/He u, O rI]/He u, and 
C iv/Hen (e.g., Hirschmann et al. 2019; Hirschmann et al. 
2022; Feltre et al. 2016). Unfortunately, many of these very- 
high ionization lines are often weak and thus may not always 
have well-constrained ratios, as is the case for the object 
studied in this work. Photoionization modeling (e.g., from 
Cleri et al. 2023), along with comparisons to data across a 
broad redshift range (0 S z x 8.5) suggest that sources of 
very-highly-ionizing photons (254.42 eV, Berg et al. 2021) 


may easily be confounded with AGN by traditional BPT-style 
diagnostics. This result has been shown with well-studied 
z ~ 0 extreme metal-poor dwarf star-forming galaxies, which 
have been used as analogs to galaxies in the EoR (e.g., Berg 
et al. 2019, 2021; Olivier et al. 2022). 

Many recent studies from early JWST data have shown ex- 
actly this: high ionization, low metallicity star formation may 
produce line ratios consistent with AGN in the lower red- 
shift diagnostics (e.g., Trump et al. 2022; Brinchmann 2022; 
Katz et al. 2022; Trussler et al. 2022). From the other side, 
other work has shown known AGNSs that produce line ra- 
tios consistent with star-forming galaxies (Übler et al. 2023). 
This growing body of work suggests that the dichotomous 
classification of a galaxy as either AGN-dominated or star- 
formation-dominated is inadequate to accurately describe a 
galaxy's ionizing spectrum in the early Universe. 

These works suggest a larger “composite” region of BPT- 
style diagrams as a function of redshift. This allows for 
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contributions to the ionizing spectrum from multiple sources 
(e.g., an accreting black hole and extreme metal-poor star 
formation). As more exotic systems in the early Universe 
are found, evolved versions of these diagnostics and their use 
with other information about the galaxy (e.g., the analysis 
of broad Balmer lines and SED fits as shown in this work) 
will become increasingly necessary to discriminate between 
sources of ionization. 


7.4. Implications for the Re-ionization of this Overdense 
Region at z = 8.7 


The presence of an AGN in this object adds another layer 
of intrigue to this region. This object is one of two spec- 
troscopically confirmed galaxies discussed by Larson et al. 
(2022), which reside within a larger photometric overdensity 
of five bright HST selected galaxies discussed by Finkelstein 
et al. (2022a). The CEERS NIRCam data also show further 
evidence of an over-abundance of z ~ 8.5-9 galaxies in this 
field, as shown in Finkelstein et al. (2022b) and studied in 
more detail by Whitler et al. (2023). Larson et al. (2022) dis- 
cussed the ability of the galaxies in this region to ionize their 
surroundings, allowing Lya to be visible. They found that 
an overdense region could provide the needed emissivity to 
reionize a large (21 pMpc) bubble, necessary for Lya to red- 
shift out of resonance. While at the observed epoch, the AGN 
does not appear to dominate the rest-UV emission, potential 
past periods of super-Eddington growth could have emitted 
large amounts of ionizing photons, potentially ionizing a large 
volume around this region. 


8. CONCLUSIONS 


We present the discovery of an accreting supermassive 
black hole at z = 8.679, using spectroscopy from NIRSpec 
and NIRCam/WFSS and imaging from NIRCam and MIRI 
from the JWST CEERS survey (Bagley et al. 2022; Finkelstein 
et al. 2022b, Finkelstein et al., prep). This source, denoted 
here as CEERS_ 1019, was initially identified asa z ~ 8 photo- 
metric Lya-break dropout candidate by Roberts-Borsani et al. 
(2016b) with spectroscopic confirmation via Lyo emission 
using Keck/MOSFIRE by Zitrin et al. (2015). We detect sev- 
eral strong rest-UV and rest-optical emission lines using the 
medium resolution JWST/NIRSpec spectroscopy (R~1000) 
covering 1—5,m. From this work, we measure a significant 
broad component of the Hf emission line with FWHM ~ 1200 
km s^!, which we conclude originates from AGN activity. 

Our measurements are based on observations from several 
JWST instruments, but our key results are derived from the 1— 
5um medium-resolution grating spectra from NIRSpec. We 
use an automated line-finding code (based on Larson et al. 
2018) to identify significant emission features in an unbiased 
and systematic way, running line-injection simulations to es- 
timate accurate line flux uncertainties. Due to the very high 


signal-to-noise of the [O m] 5008 À emission line, most of 
our emission line fits are tied to the redshift and FWHM of 
this line. We perform more customized fits when required, 
including for known doublets, as well as observed broad lines. 
We also note that a few NIRSpec-detected lines are also ob- 
servable in the NIRCam WFSS spectra, albeit at lower signal- 
to-noise. 

Our key result is that the observed HG emission line has a 
clear broad component, comprising ~ half of the emission line 
flux, with a FWHM ~ 1200 km s^!. As this broad component 
is not seen in the stronger [O m1] lines (as would be the case 
in large-scale outflows), we conclude that its origin is from a 
broad-line region around an AGN. This is supported by weak 
N V, N IV and Mg II emission, as well as weak broad C III] 
emission, and the morphology, which shows a compact point 
source amongst three Sersic-like clumps. 

We explore the properties of both this accreting SMBH 
as well as the host galaxy. Constraints from the continuum 
SED, including photometry from HST, as well as NIRCam 
and MIRI show that the continuum emission is dominated by 
stellar light, particularly in the rest-UV, and that the stellar 
population is modestly massive (log M/Mo ~ 9.5) and heavily 
star-forming (log sSFR ~ —7.9). 

From the width and flux of the broad Hf emission feature, 
we estimate the mass of the SMBH to be log (M/Mo) =6.95 
+ 0.37, and that it is accreting at 1.2 (+ 0.5) x the Eddington 
limit. From the ratios of narrow emission lines, we find that 
the gas in this galaxy is modestly metal poor (~0.1 Zo) with 
little dust attenuation, dense, and highly ionized. Similar to 
other recent JWST results, CEERS_1019 sits elevated over 
most star-forming models in diagnostic line-ratio diagrams. 
While this could be interpreted as evidence that the AGN 
significantly contributes to the narrow lines, the presence 
of other presumably star-forming dominated galaxies in a 
similar line-ratio regime means that we cannot rule out stellar- 
dominated ionization. 

We discuss the implications for the presence of this black 
hole early in cosmic history. We find that it is difficult to 
explain a SMBH of this mass at z ~ 8.7 with a stellar seed, 
unless periodic episodes of super-Eddington accretion are 
possible. Alternatively, Eddington-limited accretion from a 
massive (^10^-6 Mo) direct-collapse black hole seed could 
reach the target mass by the observed epoch. Either scenario 
is somewhat exotic — the uncovering of a larger population of 
early SMBHs will place further constraints on their seeding 
and growth mechanisms. 

We conclude by noting that while the broad H8 component 
is statistically significant (2.50) and is significantly required 
by the fit (Agrc ~ 3), we can only detect this component in 
one Balmer line. However, the next few months will see a 
MIRI spectrum obtained by the MIRI GTO team (PID 1262), 
covering Ha. If, as we have concluded, a broad-line AGN is 
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present in this source, then the upcoming MIRI Ha spectra 
should show a well-detected broad line. 
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APPENDIX 


A. FAR-INFRARED AND SUBMILLIMETER CONSTRAINTS 


Recent studies have shown that a significant fraction of high-redshift AGN exhibit large rest-frame IR luminosities, with values 
ranging 10? — 10! Lo (e.g., Decarli et al. 2018). Here we report on the far-infrared and submillimeter photometric constraints on 
this source, exploiting the deep ancillary data in the field. As shown in Figure 12, a ~ 6c SCUBA-2 detection (yellow contours) 
with a flux of S850, = 1.9 €x 0.3 mJy (Zavala et al. 2017) is found 2" away from the position of CEERS. 1019 (red circle). Zavala 
et al. 2018 associated this emission as originating from a nearby Zphot ¥ 3 galaxy (named 850.44) marked by a blue cross in Figure 
12. This nearby source is also detected with Spitzer at 24 um and at 100m with Herschel, as shown in panels 2 and 3 of Figure 
12. A recent survey was conducted with the JVLA at 3 GHz (PI: M. Dickinson; see also Jimenez-Andrade et al. in prep.), which 
also found detectable emission at this location, likely associated with the nearby source (right panel in Figure 12). 


z=8.679 AGN witH CEERS 25 


If this submm detection was instead emitting from CEERS_1019 at z = 8.679, it would imply an IR luminosity of Lip ~ 
3x 10!? Lo. Although this is in line with the luminosity of z ~ 7 quasars detected with ALMA (e.g. Decarli et al. 2018), it would 
exceed the Eddington limit of our target by a factor of ~ 10. This supports the assumption that the submillimeter emission arises 
(mainly) from a different neighboring source. 


eC F277W a MIPS 24um = PACS 1000m EA — 


AT a 


Figure 12. From left to right: JWST/NIRCam F277W, Spitzer/MIPS 24 um, Herschel/PACS 100 um, and VLA 3 GHz (north is up, east is to 
the left). All the 20" x 20" cutouts are centered on the coordinates of the 60 SCUBA-2/850 um detection (yellow contours) found around the 
location of CEERS 1019 (red circle). The blue cross marks the position of a z ~ 3 galaxy detected at 24 um, and 3 GHz, which is likely the 
correct counterpart and the main contributor to the SCUBA-2 detection. 


B. X-RAY CONSTRAINTS 


The Chandra X-ray Observatory took an 800 ks exposure over the EGS field (Nandra et al. 2015), but there is no emission 
detected at the location of CEERS_1019 in this data. Adopting a 0.5-10 keV sensitivity of 8.22 x 1076 erg cm! s^! (Nandra 
et al. 2015) and assuming a photon index (I) of 1.4, we estimate an upper limit of Lx < 107^? erg s~!. This constraint places 
CEERS 1019 around or below the knee luminosity (L5.) of the AGN X-ray luminosity function (i.e., the Seyfert regime) at lower 
redshifts (z ~ 0-5; e.g., Aird et al. 2015). Thanks to its unprecedented sensitivity, JWST begins to uncover the Seyfert-like AGN 
population in the early universe. 


Figure 13. Full-band (0.5-7 keV) X-ray images from Nandra et al. (2015). The format is similar to Fig. 12. The left panel shows the original 
image, where each white pixel indicates one (or multiple) X-ray photon(s). The right panel is a smoothed version of the left. Noise dominates 
in the region around the source. 
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